Re: [Gsll-devel] practical differences between native CL array and the GSLL marray type?
Thanks for adding this so quick Thomas! Just pulled down the new version and it looks good. The specialisation for gsll-marrays looks like it should be quite easy, I might even tackle it myself if Liam has no objections..
One thing that would also be useful to me in my current work is to compare matrices for some degree of approximate equality - eg if Y is the pseudoinverse of the matrix X, then the properties XYX = X and YXY = Y should hold, but there is likely to be some numerical aberrations. I thought about a function for numerical-approximately-equal, with a tolerance parameter, but now I think probably the best way to handle this is to pass in a custom :test function. Would you agree?
Malcolm
On Wed, Apr 15, 2009 at 5:02 AM, Thomas M. Hermann tmh.public@gmail.com wrote:
The NUMERICAL-EQUAL function has been implemented as a generic function in LISP-UNIT. It is specialized for the system classes NUMBER, LIST, VECTOR and ARRAY. It has been lightly tested. The documentation has not been updated. Please note that when you specialize it for your data type, a congruent argument list is (RESULT1 RESULT2 &KEY TEST). This will be clear when the documentation is updated.
Also note that the versions specialized for the CL system classes are for RESULT1 and RESULT2 of the same class, so, for example
(NUMERICAL-EQUAL LIST ARRAY) ; Not defined => an error
will throw an error. I'll get around to defining the mixed class versions later. Sooner if there is a real need for them.
Please feel free to email me bugs and further suggestions.
Cheers,
Tom
Liam Healy wrote:
I haven't defined such a predicate, but I think it's a good idea. As you have seen, in order to do tests on marrays, I turn all the marrays into CL arrays and do the comparisons on the result. I don't know if it's the best possible way to do it, but it works. You get an error from assert-numerical-equal because it uses #'numerical-equal as a test, and that does not handle marrays. Your suggestion spurred me to take a look at the lisp-unit source; what I'm thinking about is changing #'numerical-equal to a generic function with some methods pre-defined for CL classes; this would permit applications like GSLL to add methods for their own objects.
I've cc'ed Tom Hermann on this email; he is behind the modified version of lisp-unit more so than I am. Tom: what do you think of this idea?
Liam
On Mon, Apr 13, 2009 at 11:05 AM, Malcolm Reynolds malcolm.reynolds@gmail.com wrote:
On a related point, I'm writing some unit tests using your modified version of lisp-unit (the one from repo.or.cz) and I haven't as yet been able to figure out a good way to compare two matrices for equality. In the automatically generated tests in the GSLL code it appears you do (cl-array ...) to convert the matrices back to native form. Is this the best way possible, or is there some numerical equality predicate in the gsll library that I've missed? It would seem to be ideal to be able to define a test by (assert-numerical-equal <marray 1> <marray 2>) but I get an error message when I try this...
Malcolm
On Mon, Apr 13, 2009 at 2:24 PM, Malcolm Reynolds malcolm.reynolds@gmail.com wrote:
Hi Liam, thanks for your response.
Glad to hear that marrays are recommended, and especially that they are implemented well on SBCL as that is what I'm using.
As for the million element square matrices then yes, clearly on reflection that won't work. Some of the matrices I'm using will be sparse, but some will not be sparse at all - without wanting to get bogged down in details, I'm going to need to represent the Graph Laplacian, which is a matrix that is equal to the degree matrix minus the adjacency matrix (and so in a fairly sparse graph this matrix is sparse) but I will also need the graph Kernel, which is the pseudoinverse of the Laplacian. This will not be sparse at all, so clearly a direct representation of it isn't possible. I know of someone in my university's department who has a fast method for this kind of situation which means that calculating the full kernel matrix is unnecessary, so I guess I can look into that if and when I find that performance is a limiting factor.
Anyway thanks again for your response, and for your hard work on GSLL. I'm pretty new to lisp and I'm continuously impressed by the quality of the libraries!
Malcolm
On Sat, Apr 11, 2009 at 3:46 PM, Liam Healy lhealy@common-lisp.net wrote:
Yes, it is definitely my recommendation to use marrays everywhere. I am starting to get in the habit of doing so myself. One of my ideas is to extend CL array operations to marrays so that you need not convert back and forth. It is also my intent that other math software libraries would also use marrays (or an extension of them), so that they become a standard way to interchange data between libraries.
There is more overhead for marrays than plain arrays, but on SBCL where the native arrays are used in GSL directly, there is very little, and there is no copying. I recommend that you build it with marray and then see if you have performance problems. However, if you are talking matrices of 1000000x1000000, I think you will have problems in any representation on almost any kind of hardware. For single floats, one such matrix is 4TB. In that case, I recommend you rethink the problem structure. If these are sparse matrices, you should use a sparse matrix representation. GSL doesn't have such a representation, but others have been interested in having it so it might be worth posting something to the GSL mailing list -- perhaps someone has written a GSL-compatible sparse array representation.
Liam
On Sat, Apr 11, 2009 at 8:12 AM, Malcolm Reynolds malcolm.reynolds@gmail.com wrote:
Hi all, hope this is the right place for this kind of request.
I'm writing a program which will be doing Machine Learning over graph data structures. In it I'm going to be needing to represent a bunch of matrices, some sparse, some definitely not sparse, and do fairly basic operations on these matrices (add/ subtract/multiply, pseudoinverse, probably some others). I've started a basic prototype using the standard CL arrays but I always knew I would need GSLL to calculate a pseudo inverse.
My question is, would it be a reasonably sensible decision to just use marrays for the whole program rather than worrying about converting to/ from CL-type arrays all the time? Are there any places where CL-arrays beat marrays either in terms of memory usage, access speed, or anything else? I might be needing to scale this system up be dealing with N*N square matrices where N is on the order of hundreds of thousands or even millions, if that makes any difference to the recommendation.
Thanks in advance for any help!
Malcolm Reynolds
-- === Thomas M. Hermann ===
I certainly have no objections to any contributions. Go ahead and write a method and post it here. I would make a simple (numerical-equal (cl-array a) (cl-array b) ...) kind of test.
Liam
On Wed, Apr 15, 2009 at 10:46 AM, Malcolm Reynolds malcolm.reynolds@gmail.com wrote:
Thanks for adding this so quick Thomas! Just pulled down the new version and it looks good. The specialisation for gsll-marrays looks like it should be quite easy, I might even tackle it myself if Liam has no objections..
One thing that would also be useful to me in my current work is to compare matrices for some degree of approximate equality - eg if Y is the pseudoinverse of the matrix X, then the properties XYX = X and YXY = Y should hold, but there is likely to be some numerical aberrations. I thought about a function for numerical-approximately-equal, with a tolerance parameter, but now I think probably the best way to handle this is to pass in a custom :test function. Would you agree?
Malcolm
On Wed, Apr 15, 2009 at 5:02 AM, Thomas M. Hermann tmh.public@gmail.com wrote:
The NUMERICAL-EQUAL function has been implemented as a generic function in LISP-UNIT. It is specialized for the system classes NUMBER, LIST, VECTOR and ARRAY. It has been lightly tested. The documentation has not been updated. Please note that when you specialize it for your data type, a congruent argument list is (RESULT1 RESULT2 &KEY TEST). This will be clear when the documentation is updated.
Also note that the versions specialized for the CL system classes are for RESULT1 and RESULT2 of the same class, so, for example
(NUMERICAL-EQUAL LIST ARRAY) ; Not defined => an error
will throw an error. I'll get around to defining the mixed class versions later. Sooner if there is a real need for them.
Please feel free to email me bugs and further suggestions.
Cheers,
Tom
Liam Healy wrote:
I haven't defined such a predicate, but I think it's a good idea. As you have seen, in order to do tests on marrays, I turn all the marrays into CL arrays and do the comparisons on the result. I don't know if it's the best possible way to do it, but it works. You get an error from assert-numerical-equal because it uses #'numerical-equal as a test, and that does not handle marrays. Your suggestion spurred me to take a look at the lisp-unit source; what I'm thinking about is changing #'numerical-equal to a generic function with some methods pre-defined for CL classes; this would permit applications like GSLL to add methods for their own objects.
I've cc'ed Tom Hermann on this email; he is behind the modified version of lisp-unit more so than I am. Tom: what do you think of this idea?
Liam
On Mon, Apr 13, 2009 at 11:05 AM, Malcolm Reynolds malcolm.reynolds@gmail.com wrote:
On a related point, I'm writing some unit tests using your modified version of lisp-unit (the one from repo.or.cz) and I haven't as yet been able to figure out a good way to compare two matrices for equality. In the automatically generated tests in the GSLL code it appears you do (cl-array ...) to convert the matrices back to native form. Is this the best way possible, or is there some numerical equality predicate in the gsll library that I've missed? It would seem to be ideal to be able to define a test by (assert-numerical-equal <marray 1> <marray 2>) but I get an error message when I try this...
Malcolm
On Mon, Apr 13, 2009 at 2:24 PM, Malcolm Reynolds malcolm.reynolds@gmail.com wrote:
Hi Liam, thanks for your response.
Glad to hear that marrays are recommended, and especially that they are implemented well on SBCL as that is what I'm using.
As for the million element square matrices then yes, clearly on reflection that won't work. Some of the matrices I'm using will be sparse, but some will not be sparse at all - without wanting to get bogged down in details, I'm going to need to represent the Graph Laplacian, which is a matrix that is equal to the degree matrix minus the adjacency matrix (and so in a fairly sparse graph this matrix is sparse) but I will also need the graph Kernel, which is the pseudoinverse of the Laplacian. This will not be sparse at all, so clearly a direct representation of it isn't possible. I know of someone in my university's department who has a fast method for this kind of situation which means that calculating the full kernel matrix is unnecessary, so I guess I can look into that if and when I find that performance is a limiting factor.
Anyway thanks again for your response, and for your hard work on GSLL. I'm pretty new to lisp and I'm continuously impressed by the quality of the libraries!
Malcolm
On Sat, Apr 11, 2009 at 3:46 PM, Liam Healy lhealy@common-lisp.net wrote:
Yes, it is definitely my recommendation to use marrays everywhere. I am starting to get in the habit of doing so myself. One of my ideas is to extend CL array operations to marrays so that you need not convert back and forth. It is also my intent that other math software libraries would also use marrays (or an extension of them), so that they become a standard way to interchange data between libraries.
There is more overhead for marrays than plain arrays, but on SBCL where the native arrays are used in GSL directly, there is very little, and there is no copying. I recommend that you build it with marray and then see if you have performance problems. However, if you are talking matrices of 1000000x1000000, I think you will have problems in any representation on almost any kind of hardware. For single floats, one such matrix is 4TB. In that case, I recommend you rethink the problem structure. If these are sparse matrices, you should use a sparse matrix representation. GSL doesn't have such a representation, but others have been interested in having it so it might be worth posting something to the GSL mailing list -- perhaps someone has written a GSL-compatible sparse array representation.
Liam
On Sat, Apr 11, 2009 at 8:12 AM, Malcolm Reynolds malcolm.reynolds@gmail.com wrote: > Hi all, hope this is the right place for this kind of request. > > I'm writing a program which will be doing Machine Learning over graph > data structures. In it I'm going to be needing to represent a bunch of > matrices, some sparse, some definitely not sparse, and do fairly basic > operations on these matrices (add/ subtract/multiply, pseudoinverse, > probably some others). I've started a basic prototype using the > standard CL arrays but I always knew I would need GSLL to calculate a > pseudo inverse. > > My question is, would it be a reasonably sensible decision to just use > marrays for the whole program rather than worrying > about converting to/ from CL-type arrays all the time? Are there any > places where CL-arrays beat marrays either in terms of memory usage, > access speed, or anything else? I might be needing to scale this > system up be dealing with N*N square matrices where N is on the order > of hundreds of thousands or even millions, if that makes any > difference to the recommendation. > > Thanks in advance for any help! > > Malcolm Reynolds
-- === Thomas M. Hermann ===
This is what I put together quickly, seemed to work for the few cases I tried it on. It will return nil when you give it a matrix and a vector that are equal in all their values, ie the vector '(1 2 3) and the matrix '((1 2 3)) won't be equal.. I'm not sure whether that's a good or bad behaviour at this point...
(defmethod lisp-unit:numerical-equal ((result1 marray) (result2 marray) &key (test #'lisp-unit:number-equal)) "Return true if the arrays are numerically equal according to :TEST." (when (equal (dimensions result1) (dimensions result2)) (let ((m (dim0 result1)) (n (dim1 result1))) (do ((i 0 (1+ i))) ((= i m)) (do ((j 0 (1+ j))) ((= j n)) (when (not (funcall test (maref result1 i j) (maref result2 i j))) (return-from lisp-unit:numerical-equal nil)))) t)))
It could definitely be a lot more pretty, with something like the do-matrix macro I pasted in #lisp, but from the comments from tcr it seems like that macro has some subtle flaws which mean it definitely isn't ready to go in any libraries for general consumption yet.
On Wed, Apr 15, 2009 at 11:36 PM, Liam Healy lhealy@common-lisp.net wrote:
I certainly have no objections to any contributions. Go ahead and write a method and post it here. I would make a simple (numerical-equal (cl-array a) (cl-array b) ...) kind of test.
Liam
On Wed, Apr 15, 2009 at 10:46 AM, Malcolm Reynolds malcolm.reynolds@gmail.com wrote:
Thanks for adding this so quick Thomas! Just pulled down the new version and it looks good. The specialisation for gsll-marrays looks like it should be quite easy, I might even tackle it myself if Liam has no objections..
One thing that would also be useful to me in my current work is to compare matrices for some degree of approximate equality - eg if Y is the pseudoinverse of the matrix X, then the properties XYX = X and YXY = Y should hold, but there is likely to be some numerical aberrations. I thought about a function for numerical-approximately-equal, with a tolerance parameter, but now I think probably the best way to handle this is to pass in a custom :test function. Would you agree?
Malcolm
On Wed, Apr 15, 2009 at 5:02 AM, Thomas M. Hermann tmh.public@gmail.com wrote:
The NUMERICAL-EQUAL function has been implemented as a generic function in LISP-UNIT. It is specialized for the system classes NUMBER, LIST, VECTOR and ARRAY. It has been lightly tested. The documentation has not been updated. Please note that when you specialize it for your data type, a congruent argument list is (RESULT1 RESULT2 &KEY TEST). This will be clear when the documentation is updated.
Also note that the versions specialized for the CL system classes are for RESULT1 and RESULT2 of the same class, so, for example
(NUMERICAL-EQUAL LIST ARRAY) ; Not defined => an error
will throw an error. I'll get around to defining the mixed class versions later. Sooner if there is a real need for them.
Please feel free to email me bugs and further suggestions.
Cheers,
Tom
Liam Healy wrote:
I haven't defined such a predicate, but I think it's a good idea. As you have seen, in order to do tests on marrays, I turn all the marrays into CL arrays and do the comparisons on the result. I don't know if it's the best possible way to do it, but it works. You get an error from assert-numerical-equal because it uses #'numerical-equal as a test, and that does not handle marrays. Your suggestion spurred me to take a look at the lisp-unit source; what I'm thinking about is changing #'numerical-equal to a generic function with some methods pre-defined for CL classes; this would permit applications like GSLL to add methods for their own objects.
I've cc'ed Tom Hermann on this email; he is behind the modified version of lisp-unit more so than I am. Tom: what do you think of this idea?
Liam
On Mon, Apr 13, 2009 at 11:05 AM, Malcolm Reynolds malcolm.reynolds@gmail.com wrote:
On a related point, I'm writing some unit tests using your modified version of lisp-unit (the one from repo.or.cz) and I haven't as yet been able to figure out a good way to compare two matrices for equality. In the automatically generated tests in the GSLL code it appears you do (cl-array ...) to convert the matrices back to native form. Is this the best way possible, or is there some numerical equality predicate in the gsll library that I've missed? It would seem to be ideal to be able to define a test by (assert-numerical-equal <marray 1> <marray 2>) but I get an error message when I try this...
Malcolm
On Mon, Apr 13, 2009 at 2:24 PM, Malcolm Reynolds malcolm.reynolds@gmail.com wrote:
Hi Liam, thanks for your response.
Glad to hear that marrays are recommended, and especially that they are implemented well on SBCL as that is what I'm using.
As for the million element square matrices then yes, clearly on reflection that won't work. Some of the matrices I'm using will be sparse, but some will not be sparse at all - without wanting to get bogged down in details, I'm going to need to represent the Graph Laplacian, which is a matrix that is equal to the degree matrix minus the adjacency matrix (and so in a fairly sparse graph this matrix is sparse) but I will also need the graph Kernel, which is the pseudoinverse of the Laplacian. This will not be sparse at all, so clearly a direct representation of it isn't possible. I know of someone in my university's department who has a fast method for this kind of situation which means that calculating the full kernel matrix is unnecessary, so I guess I can look into that if and when I find that performance is a limiting factor.
Anyway thanks again for your response, and for your hard work on GSLL. I'm pretty new to lisp and I'm continuously impressed by the quality of the libraries!
Malcolm
On Sat, Apr 11, 2009 at 3:46 PM, Liam Healy lhealy@common-lisp.net wrote: > Yes, it is definitely my recommendation to use marrays > everywhere. I am starting to get in the habit of doing so > myself. One of my ideas is to extend CL array operations to > marrays so that you need not convert back and forth. It is > also my intent that other math software libraries would also > use marrays (or an extension of them), so that they become a > standard way to interchange data between libraries. > > There is more overhead for marrays than plain arrays, but on > SBCL where the native arrays are used in GSL directly, there > is very little, and there is no copying. I recommend that > you build it with marray and then see if you have > performance problems. However, if you are talking matrices > of 1000000x1000000, I think you will have problems in any > representation on almost any kind of hardware. For single > floats, one such matrix is 4TB. In that case, I recommend > you rethink the problem structure. If these are sparse > matrices, you should use a sparse matrix representation. > GSL doesn't have such a representation, but others have been > interested in having it so it might be worth posting > something to the GSL mailing list -- perhaps someone has > written a GSL-compatible sparse array representation. > > Liam > > > > On Sat, Apr 11, 2009 at 8:12 AM, Malcolm Reynolds > malcolm.reynolds@gmail.com wrote: >> Hi all, hope this is the right place for this kind of request. >> >> I'm writing a program which will be doing Machine Learning over graph >> data structures. In it I'm going to be needing to represent a bunch of >> matrices, some sparse, some definitely not sparse, and do fairly basic >> operations on these matrices (add/ subtract/multiply, pseudoinverse, >> probably some others). I've started a basic prototype using the >> standard CL arrays but I always knew I would need GSLL to calculate a >> pseudo inverse. >> >> My question is, would it be a reasonably sensible decision to just use >> marrays for the whole program rather than worrying >> about converting to/ from CL-type arrays all the time? Are there any >> places where CL-arrays beat marrays either in terms of memory usage, >> access speed, or anything else? I might be needing to scale this >> system up be dealing with N*N square matrices where N is on the order >> of hundreds of thousands or even millions, if that makes any >> difference to the recommendation. >> >> Thanks in advance for any help! >> >> Malcolm Reynolds
-- === Thomas M. Hermann ===
I'd make it simpler:
(defmethod lisp-unit:numerical-equal ((result1 marray) (result2 marray) &key (test #'lisp-unit:number-equal)) "Return true if the arrays are numerically equal according to :TEST." (when (equal (dimensions result1) (dimensions result2)) (lisp-unit:numerical-equal (cl-array result1) (cl-array result2)) :test test))
works?
Liam
On Wed, Apr 15, 2009 at 7:06 PM, Malcolm Reynolds malcolm.reynolds@gmail.com wrote:
This is what I put together quickly, seemed to work for the few cases I tried it on. It will return nil when you give it a matrix and a vector that are equal in all their values, ie the vector '(1 2 3) and the matrix '((1 2 3)) won't be equal.. I'm not sure whether that's a good or bad behaviour at this point...
(defmethod lisp-unit:numerical-equal ((result1 marray) (result2 marray) &key (test #'lisp-unit:number-equal)) "Return true if the arrays are numerically equal according to :TEST." (when (equal (dimensions result1) (dimensions result2)) (let ((m (dim0 result1)) (n (dim1 result1))) (do ((i 0 (1+ i))) ((= i m)) (do ((j 0 (1+ j))) ((= j n)) (when (not (funcall test (maref result1 i j) (maref result2 i j))) (return-from lisp-unit:numerical-equal nil)))) t)))
It could definitely be a lot more pretty, with something like the do-matrix macro I pasted in #lisp, but from the comments from tcr it seems like that macro has some subtle flaws which mean it definitely isn't ready to go in any libraries for general consumption yet.
On Wed, Apr 15, 2009 at 11:36 PM, Liam Healy lhealy@common-lisp.net wrote:
I certainly have no objections to any contributions. Go ahead and write a method and post it here. I would make a simple (numerical-equal (cl-array a) (cl-array b) ...) kind of test.
Liam
On Wed, Apr 15, 2009 at 10:46 AM, Malcolm Reynolds malcolm.reynolds@gmail.com wrote:
Thanks for adding this so quick Thomas! Just pulled down the new version and it looks good. The specialisation for gsll-marrays looks like it should be quite easy, I might even tackle it myself if Liam has no objections..
One thing that would also be useful to me in my current work is to compare matrices for some degree of approximate equality - eg if Y is the pseudoinverse of the matrix X, then the properties XYX = X and YXY = Y should hold, but there is likely to be some numerical aberrations. I thought about a function for numerical-approximately-equal, with a tolerance parameter, but now I think probably the best way to handle this is to pass in a custom :test function. Would you agree?
Malcolm
On Wed, Apr 15, 2009 at 5:02 AM, Thomas M. Hermann tmh.public@gmail.com wrote:
The NUMERICAL-EQUAL function has been implemented as a generic function in LISP-UNIT. It is specialized for the system classes NUMBER, LIST, VECTOR and ARRAY. It has been lightly tested. The documentation has not been updated. Please note that when you specialize it for your data type, a congruent argument list is (RESULT1 RESULT2 &KEY TEST). This will be clear when the documentation is updated.
Also note that the versions specialized for the CL system classes are for RESULT1 and RESULT2 of the same class, so, for example
(NUMERICAL-EQUAL LIST ARRAY) ; Not defined => an error
will throw an error. I'll get around to defining the mixed class versions later. Sooner if there is a real need for them.
Please feel free to email me bugs and further suggestions.
Cheers,
Tom
Liam Healy wrote:
I haven't defined such a predicate, but I think it's a good idea. As you have seen, in order to do tests on marrays, I turn all the marrays into CL arrays and do the comparisons on the result. I don't know if it's the best possible way to do it, but it works. You get an error from assert-numerical-equal because it uses #'numerical-equal as a test, and that does not handle marrays. Your suggestion spurred me to take a look at the lisp-unit source; what I'm thinking about is changing #'numerical-equal to a generic function with some methods pre-defined for CL classes; this would permit applications like GSLL to add methods for their own objects.
I've cc'ed Tom Hermann on this email; he is behind the modified version of lisp-unit more so than I am. Tom: what do you think of this idea?
Liam
On Mon, Apr 13, 2009 at 11:05 AM, Malcolm Reynolds malcolm.reynolds@gmail.com wrote:
On a related point, I'm writing some unit tests using your modified version of lisp-unit (the one from repo.or.cz) and I haven't as yet been able to figure out a good way to compare two matrices for equality. In the automatically generated tests in the GSLL code it appears you do (cl-array ...) to convert the matrices back to native form. Is this the best way possible, or is there some numerical equality predicate in the gsll library that I've missed? It would seem to be ideal to be able to define a test by (assert-numerical-equal <marray 1> <marray 2>) but I get an error message when I try this...
Malcolm
On Mon, Apr 13, 2009 at 2:24 PM, Malcolm Reynolds malcolm.reynolds@gmail.com wrote: > Hi Liam, thanks for your response. > > Glad to hear that marrays are recommended, and especially that they > are implemented well on SBCL as that is what I'm using. > > As for the million element square matrices then yes, clearly on > reflection that won't work. Some of the matrices I'm using will be > sparse, but some will not be sparse at all - without wanting to get > bogged down in details, I'm going to need to represent the Graph > Laplacian, which is a matrix that is equal to the degree matrix minus > the adjacency matrix (and so in a fairly sparse graph this matrix is > sparse) but I will also need the graph Kernel, which is the > pseudoinverse of the Laplacian. This will not be sparse at all, so > clearly a direct representation of it isn't possible. I know of > someone in my university's department who has a fast method for this > kind of situation which means that calculating the full kernel matrix > is unnecessary, so I guess I can look into that if and when I find > that performance is a limiting factor. > > Anyway thanks again for your response, and for your hard work on GSLL. > I'm pretty new to lisp and I'm continuously impressed by the quality > of the libraries! > > Malcolm > > On Sat, Apr 11, 2009 at 3:46 PM, Liam Healy lhealy@common-lisp.net wrote: >> Yes, it is definitely my recommendation to use marrays >> everywhere. I am starting to get in the habit of doing so >> myself. One of my ideas is to extend CL array operations to >> marrays so that you need not convert back and forth. It is >> also my intent that other math software libraries would also >> use marrays (or an extension of them), so that they become a >> standard way to interchange data between libraries. >> >> There is more overhead for marrays than plain arrays, but on >> SBCL where the native arrays are used in GSL directly, there >> is very little, and there is no copying. I recommend that >> you build it with marray and then see if you have >> performance problems. However, if you are talking matrices >> of 1000000x1000000, I think you will have problems in any >> representation on almost any kind of hardware. For single >> floats, one such matrix is 4TB. In that case, I recommend >> you rethink the problem structure. If these are sparse >> matrices, you should use a sparse matrix representation. >> GSL doesn't have such a representation, but others have been >> interested in having it so it might be worth posting >> something to the GSL mailing list -- perhaps someone has >> written a GSL-compatible sparse array representation. >> >> Liam >> >> >> >> On Sat, Apr 11, 2009 at 8:12 AM, Malcolm Reynolds >> malcolm.reynolds@gmail.com wrote: >>> Hi all, hope this is the right place for this kind of request. >>> >>> I'm writing a program which will be doing Machine Learning over graph >>> data structures. In it I'm going to be needing to represent a bunch of >>> matrices, some sparse, some definitely not sparse, and do fairly basic >>> operations on these matrices (add/ subtract/multiply, pseudoinverse, >>> probably some others). I've started a basic prototype using the >>> standard CL arrays but I always knew I would need GSLL to calculate a >>> pseudo inverse. >>> >>> My question is, would it be a reasonably sensible decision to just use >>> marrays for the whole program rather than worrying >>> about converting to/ from CL-type arrays all the time? Are there any >>> places where CL-arrays beat marrays either in terms of memory usage, >>> access speed, or anything else? I might be needing to scale this >>> system up be dealing with N*N square matrices where N is on the order >>> of hundreds of thousands or even millions, if that makes any >>> difference to the recommendation. >>> >>> Thanks in advance for any help! >>> >>> Malcolm Reynolds
-- === Thomas M. Hermann ===
Yeah, I guess that is fine too but I was attempting to avoid using (cl-array ..). I know you said on SBCL there is low overhead for the conversion but I figured in case that didn't hold on other CL implementations this probably does. Since it's now a generic function it seemed worth making it specific on your custom datatypes, but I appreciate what seems like it should be more optimal might not always turn out to be so..
Malcolm
On Thu, Apr 16, 2009 at 12:28 AM, Liam Healy lhealy@common-lisp.net wrote:
I'd make it simpler:
(defmethod lisp-unit:numerical-equal ((result1 marray) (result2 marray) &key (test #'lisp-unit:number-equal)) "Return true if the arrays are numerically equal according to :TEST." (when (equal (dimensions result1) (dimensions result2)) (lisp-unit:numerical-equal (cl-array result1) (cl-array result2)) :test test))
works?
Liam
On Wed, Apr 15, 2009 at 7:06 PM, Malcolm Reynolds malcolm.reynolds@gmail.com wrote:
This is what I put together quickly, seemed to work for the few cases I tried it on. It will return nil when you give it a matrix and a vector that are equal in all their values, ie the vector '(1 2 3) and the matrix '((1 2 3)) won't be equal.. I'm not sure whether that's a good or bad behaviour at this point...
(defmethod lisp-unit:numerical-equal ((result1 marray) (result2 marray) &key (test #'lisp-unit:number-equal)) "Return true if the arrays are numerically equal according to :TEST." (when (equal (dimensions result1) (dimensions result2)) (let ((m (dim0 result1)) (n (dim1 result1))) (do ((i 0 (1+ i))) ((= i m)) (do ((j 0 (1+ j))) ((= j n)) (when (not (funcall test (maref result1 i j) (maref result2 i j))) (return-from lisp-unit:numerical-equal nil)))) t)))
It could definitely be a lot more pretty, with something like the do-matrix macro I pasted in #lisp, but from the comments from tcr it seems like that macro has some subtle flaws which mean it definitely isn't ready to go in any libraries for general consumption yet.
On Wed, Apr 15, 2009 at 11:36 PM, Liam Healy lhealy@common-lisp.net wrote:
I certainly have no objections to any contributions. Go ahead and write a method and post it here. I would make a simple (numerical-equal (cl-array a) (cl-array b) ...) kind of test.
Liam
On Wed, Apr 15, 2009 at 10:46 AM, Malcolm Reynolds malcolm.reynolds@gmail.com wrote:
Thanks for adding this so quick Thomas! Just pulled down the new version and it looks good. The specialisation for gsll-marrays looks like it should be quite easy, I might even tackle it myself if Liam has no objections..
One thing that would also be useful to me in my current work is to compare matrices for some degree of approximate equality - eg if Y is the pseudoinverse of the matrix X, then the properties XYX = X and YXY = Y should hold, but there is likely to be some numerical aberrations. I thought about a function for numerical-approximately-equal, with a tolerance parameter, but now I think probably the best way to handle this is to pass in a custom :test function. Would you agree?
Malcolm
On Wed, Apr 15, 2009 at 5:02 AM, Thomas M. Hermann tmh.public@gmail.com wrote:
The NUMERICAL-EQUAL function has been implemented as a generic function in LISP-UNIT. It is specialized for the system classes NUMBER, LIST, VECTOR and ARRAY. It has been lightly tested. The documentation has not been updated. Please note that when you specialize it for your data type, a congruent argument list is (RESULT1 RESULT2 &KEY TEST). This will be clear when the documentation is updated.
Also note that the versions specialized for the CL system classes are for RESULT1 and RESULT2 of the same class, so, for example
(NUMERICAL-EQUAL LIST ARRAY) ; Not defined => an error
will throw an error. I'll get around to defining the mixed class versions later. Sooner if there is a real need for them.
Please feel free to email me bugs and further suggestions.
Cheers,
Tom
Liam Healy wrote:
I haven't defined such a predicate, but I think it's a good idea. As you have seen, in order to do tests on marrays, I turn all the marrays into CL arrays and do the comparisons on the result. I don't know if it's the best possible way to do it, but it works. You get an error from assert-numerical-equal because it uses #'numerical-equal as a test, and that does not handle marrays. Your suggestion spurred me to take a look at the lisp-unit source; what I'm thinking about is changing #'numerical-equal to a generic function with some methods pre-defined for CL classes; this would permit applications like GSLL to add methods for their own objects.
I've cc'ed Tom Hermann on this email; he is behind the modified version of lisp-unit more so than I am. Tom: what do you think of this idea?
Liam
On Mon, Apr 13, 2009 at 11:05 AM, Malcolm Reynolds malcolm.reynolds@gmail.com wrote: > On a related point, I'm writing some unit tests using your modified > version of lisp-unit (the one from repo.or.cz) and I haven't as yet > been able to figure out a good way to compare two matrices for > equality. In the automatically generated tests in the GSLL code it > appears you do (cl-array ...) to convert the matrices back to native > form. Is this the best way possible, or is there some numerical > equality predicate in the gsll library that I've missed? It would seem > to be ideal to be able to > define a test by (assert-numerical-equal <marray 1> <marray 2>) but I > get an error message when I try this... > > Malcolm > > On Mon, Apr 13, 2009 at 2:24 PM, Malcolm Reynolds > malcolm.reynolds@gmail.com wrote: >> Hi Liam, thanks for your response. >> >> Glad to hear that marrays are recommended, and especially that they >> are implemented well on SBCL as that is what I'm using. >> >> As for the million element square matrices then yes, clearly on >> reflection that won't work. Some of the matrices I'm using will be >> sparse, but some will not be sparse at all - without wanting to get >> bogged down in details, I'm going to need to represent the Graph >> Laplacian, which is a matrix that is equal to the degree matrix minus >> the adjacency matrix (and so in a fairly sparse graph this matrix is >> sparse) but I will also need the graph Kernel, which is the >> pseudoinverse of the Laplacian. This will not be sparse at all, so >> clearly a direct representation of it isn't possible. I know of >> someone in my university's department who has a fast method for this >> kind of situation which means that calculating the full kernel matrix >> is unnecessary, so I guess I can look into that if and when I find >> that performance is a limiting factor. >> >> Anyway thanks again for your response, and for your hard work on GSLL. >> I'm pretty new to lisp and I'm continuously impressed by the quality >> of the libraries! >> >> Malcolm >> >> On Sat, Apr 11, 2009 at 3:46 PM, Liam Healy lhealy@common-lisp.net wrote: >>> Yes, it is definitely my recommendation to use marrays >>> everywhere. I am starting to get in the habit of doing so >>> myself. One of my ideas is to extend CL array operations to >>> marrays so that you need not convert back and forth. It is >>> also my intent that other math software libraries would also >>> use marrays (or an extension of them), so that they become a >>> standard way to interchange data between libraries. >>> >>> There is more overhead for marrays than plain arrays, but on >>> SBCL where the native arrays are used in GSL directly, there >>> is very little, and there is no copying. I recommend that >>> you build it with marray and then see if you have >>> performance problems. However, if you are talking matrices >>> of 1000000x1000000, I think you will have problems in any >>> representation on almost any kind of hardware. For single >>> floats, one such matrix is 4TB. In that case, I recommend >>> you rethink the problem structure. If these are sparse >>> matrices, you should use a sparse matrix representation. >>> GSL doesn't have such a representation, but others have been >>> interested in having it so it might be worth posting >>> something to the GSL mailing list -- perhaps someone has >>> written a GSL-compatible sparse array representation. >>> >>> Liam >>> >>> >>> >>> On Sat, Apr 11, 2009 at 8:12 AM, Malcolm Reynolds >>> malcolm.reynolds@gmail.com wrote: >>>> Hi all, hope this is the right place for this kind of request. >>>> >>>> I'm writing a program which will be doing Machine Learning over graph >>>> data structures. In it I'm going to be needing to represent a bunch of >>>> matrices, some sparse, some definitely not sparse, and do fairly basic >>>> operations on these matrices (add/ subtract/multiply, pseudoinverse, >>>> probably some others). I've started a basic prototype using the >>>> standard CL arrays but I always knew I would need GSLL to calculate a >>>> pseudo inverse. >>>> >>>> My question is, would it be a reasonably sensible decision to just use >>>> marrays for the whole program rather than worrying >>>> about converting to/ from CL-type arrays all the time? Are there any >>>> places where CL-arrays beat marrays either in terms of memory usage, >>>> access speed, or anything else? I might be needing to scale this >>>> system up be dealing with N*N square matrices where N is on the order >>>> of hundreds of thousands or even millions, if that makes any >>>> difference to the recommendation. >>>> >>>> Thanks in advance for any help! >>>> >>>> Malcolm Reynolds
-- === Thomas M. Hermann ===
There is no difference in what your form does and what mine does in terms of efficiency, on any platform. They are identical. The first time a maref is called on a non-native implementation, the marray is checked for cl-invalid. If it is T, the entire array is copied from the C side to the CL side. The exact same thing happens when cl-array is called.
On Wed, Apr 15, 2009 at 7:40 PM, Malcolm Reynolds malcolm.reynolds@gmail.com wrote:
Yeah, I guess that is fine too but I was attempting to avoid using (cl-array ..). I know you said on SBCL there is low overhead for the conversion but I figured in case that didn't hold on other CL implementations this probably does. Since it's now a generic function it seemed worth making it specific on your custom datatypes, but I appreciate what seems like it should be more optimal might not always turn out to be so..
Malcolm
On Thu, Apr 16, 2009 at 12:28 AM, Liam Healy lhealy@common-lisp.net wrote:
I'd make it simpler:
(defmethod lisp-unit:numerical-equal ((result1 marray) (result2 marray) &key (test #'lisp-unit:number-equal)) "Return true if the arrays are numerically equal according to :TEST." (when (equal (dimensions result1) (dimensions result2)) (lisp-unit:numerical-equal (cl-array result1) (cl-array result2)) :test test))
works?
Liam
On Wed, Apr 15, 2009 at 7:06 PM, Malcolm Reynolds malcolm.reynolds@gmail.com wrote:
This is what I put together quickly, seemed to work for the few cases I tried it on. It will return nil when you give it a matrix and a vector that are equal in all their values, ie the vector '(1 2 3) and the matrix '((1 2 3)) won't be equal.. I'm not sure whether that's a good or bad behaviour at this point...
(defmethod lisp-unit:numerical-equal ((result1 marray) (result2 marray) &key (test #'lisp-unit:number-equal)) "Return true if the arrays are numerically equal according to :TEST." (when (equal (dimensions result1) (dimensions result2)) (let ((m (dim0 result1)) (n (dim1 result1))) (do ((i 0 (1+ i))) ((= i m)) (do ((j 0 (1+ j))) ((= j n)) (when (not (funcall test (maref result1 i j) (maref result2 i j))) (return-from lisp-unit:numerical-equal nil)))) t)))
It could definitely be a lot more pretty, with something like the do-matrix macro I pasted in #lisp, but from the comments from tcr it seems like that macro has some subtle flaws which mean it definitely isn't ready to go in any libraries for general consumption yet.
On Wed, Apr 15, 2009 at 11:36 PM, Liam Healy lhealy@common-lisp.net wrote:
I certainly have no objections to any contributions. Go ahead and write a method and post it here. I would make a simple (numerical-equal (cl-array a) (cl-array b) ...) kind of test.
Liam
On Wed, Apr 15, 2009 at 10:46 AM, Malcolm Reynolds malcolm.reynolds@gmail.com wrote:
Thanks for adding this so quick Thomas! Just pulled down the new version and it looks good. The specialisation for gsll-marrays looks like it should be quite easy, I might even tackle it myself if Liam has no objections..
One thing that would also be useful to me in my current work is to compare matrices for some degree of approximate equality - eg if Y is the pseudoinverse of the matrix X, then the properties XYX = X and YXY = Y should hold, but there is likely to be some numerical aberrations. I thought about a function for numerical-approximately-equal, with a tolerance parameter, but now I think probably the best way to handle this is to pass in a custom :test function. Would you agree?
Malcolm
On Wed, Apr 15, 2009 at 5:02 AM, Thomas M. Hermann tmh.public@gmail.com wrote:
The NUMERICAL-EQUAL function has been implemented as a generic function in LISP-UNIT. It is specialized for the system classes NUMBER, LIST, VECTOR and ARRAY. It has been lightly tested. The documentation has not been updated. Please note that when you specialize it for your data type, a congruent argument list is (RESULT1 RESULT2 &KEY TEST). This will be clear when the documentation is updated.
Also note that the versions specialized for the CL system classes are for RESULT1 and RESULT2 of the same class, so, for example
(NUMERICAL-EQUAL LIST ARRAY) ; Not defined => an error
will throw an error. I'll get around to defining the mixed class versions later. Sooner if there is a real need for them.
Please feel free to email me bugs and further suggestions.
Cheers,
Tom
Liam Healy wrote: > I haven't defined such a predicate, but I think it's a good idea. As > you have seen, in order to do tests on marrays, I turn all the marrays > into CL arrays and do the comparisons on the result. I don't know if > it's the best possible way to do it, but it works. You get an error > from assert-numerical-equal because it uses #'numerical-equal as a > test, and that does not handle marrays. Your suggestion spurred me > to take a look at the lisp-unit source; what I'm thinking about is > changing #'numerical-equal to a generic function with some methods > pre-defined for CL classes; this would permit applications like GSLL > to add methods for their own objects. > > I've cc'ed Tom Hermann on this email; he is behind the modified > version of lisp-unit more so than I am. > Tom: what do you think of this idea? > > Liam
I have added a method for lisp-unit:numerical-equal on marrays. I changed only the vector-add test to take advantage of the ability to now compare marrays directly, to make sure it works. It actually is a fair amount of work to rewrite the tests, so I'll probably do this gradually. Thank you Malcolm and Tom for your contributions.
P.S. Further thoughts on the efficiency of fobbing the comparison off to CL behind the method: If element values were extracted using the GSL interface (like gsl_vector_get) then copying could be avoided. However, with native arrays (for SBCL), it seemed simplest to always access elements on the CL side, and just make sure they were synchronized if on a non-native implementation. It's possible that it would be faster to use the C values but since you have to read all the elements anyway to do a comparison, I don't think it would gain that much. In any case, it would be handy to have an iterate over marray elements, another thing on my to-do list.
Liam
Unless I'm mistaken I think there is a bug in the implementation currently.. this is with the latest from the lisp-unit and gsll git repositories, such that lisp-unit is at commit e7c4faa8baf9d071972a66c62671001a62f3cc1c ("Implementations of NUMERICAL-EQUAL for mixed list/vector arguments.") and gsll is at commit 64a0b6c271530b298ce97f168d03ab715fb80a39 ("Tests marrays directly"). Full log of a slime session follows:
; SLIME 2009-03-09 CL-USER> (asdf:operate 'asdf:load-op :gsll) ; loading system definition from ; /Users/malc/opt/lisp/asdf-registry/gsll.asd into #<PACKAGE "ASDF0"> ; registering #<SYSTEM "gsll" {12295D61}> as gsll ; loading system definition from ; /Users/malc/opt/lisp/sbcl/sbcl-inst-1.0.23-x86-darwin/lib/sbcl/site-systems/trivial-garbage.asd ; into #<PACKAGE "ASDF0"> ; registering #<SYSTEM TRIVIAL-GARBAGE {124BF099}> as TRIVIAL-GARBAGE ; registering #<SYSTEM TRIVIAL-GARBAGE-TESTS {1272DDF1}> as ; TRIVIAL-GARBAGE-TESTS ; loading system definition from ; /Users/malc/opt/lisp/sbcl/sbcl-inst-1.0.23-x86-darwin/lib/sbcl/site-systems/cffi.asd ; into #<PACKAGE "ASDF0"> ; registering #<SYSTEM CFFI {119164D1}> as CFFI ; loading system definition from ; /Users/malc/opt/lisp/sbcl/sbcl-inst-1.0.23-x86-darwin/lib/sbcl/site-systems/babel.asd ; into #<PACKAGE "ASDF0"> ; registering #<SYSTEM BABEL {11CA1489}> as BABEL ; loading system definition from ; /Users/malc/opt/lisp/sbcl/sbcl-inst-1.0.23-x86-darwin/lib/sbcl/site-systems/alexandria.asd ; into #<PACKAGE "ASDF0"> ; registering #<SYSTEM :ALEXANDRIA {1276FF09}> as ALEXANDRIA ; loading system definition from ; /Users/malc/opt/lisp/sbcl/sbcl-inst-1.0.23-x86-darwin/lib/sbcl/site-systems/trivial-features.asd ; into #<PACKAGE "ASDF0"> ; registering #<SYSTEM TRIVIAL-FEATURES {129ED569}> as TRIVIAL-FEATURES NIL CL-USER> (asdf:operate 'asdf:load-op :gsll-tests) ; loading system definition from ; /Users/malc/opt/lisp/asdf-registry/gsll-tests.asd into #<PACKAGE "ASDF0"> ; registering #<SYSTEM "gsll-tests" {11B966E9}> as gsll-tests ; loading system definition from ; /Users/malc/opt/lisp/asdf-registry/lisp-unit.asd into #<PACKAGE "ASDF0"> ; registering #<SYSTEM :LISP-UNIT {12484CB9}> as LISP-UNIT NIL CL-USER> (in-package :gsll) #<PACKAGE "GSLL"> GSL> (lisp-unit:assert-numerical-equal '(1 2 3) '(1 2 4)) '(1 2 4) failed: Expected (1 2 3) but saw (1 2 4) NIL GSL> (lisp-unit:assert-numerical-equal (make-marray 'double-float :initial-contents '(1 2 3)) (make-marray 'double-float :initial-contents '(1 2 4))) T
Can you guys reproduce this? I'm SBCL 1.0.23 on OS X if that makes a difference..
Cheers
Malcolm
Additionally, if I wrap the two calls to make-marray in (cl-array ...) then the test fails correctly and says
(CL-ARRAY (MAKE-MARRAY 'DOUBLE-FLOAT :INITIAL-CONTENTS '(1 2 4))) failed: Expected #(1.0d0 2.0d0 3.0d0) but saw #(1.0d0 2.0d0 4.0d0)
Malcolm
Yikes, that's an embarrassing goof in the lisp-unit:numerical-equal method for marrays. Should be fixed now, try a fresh pull and see what happens. Thanks for the report.
Liam
On Thu, Apr 30, 2009 at 10:46 AM, Malcolm Reynolds malcolm.reynolds@gmail.com wrote:
Unless I'm mistaken I think there is a bug in the implementation currently.. this is with the latest from the lisp-unit and gsll git repositories, such that lisp-unit is at commit e7c4faa8baf9d071972a66c62671001a62f3cc1c ("Implementations of NUMERICAL-EQUAL for mixed list/vector arguments.") and gsll is at commit 64a0b6c271530b298ce97f168d03ab715fb80a39 ("Tests marrays directly"). Full log of a slime session follows:
; SLIME 2009-03-09 CL-USER> (asdf:operate 'asdf:load-op :gsll) ; loading system definition from ; /Users/malc/opt/lisp/asdf-registry/gsll.asd into #<PACKAGE "ASDF0"> ; registering #<SYSTEM "gsll" {12295D61}> as gsll ; loading system definition from ; /Users/malc/opt/lisp/sbcl/sbcl-inst-1.0.23-x86-darwin/lib/sbcl/site-systems/trivial-garbage.asd ; into #<PACKAGE "ASDF0"> ; registering #<SYSTEM TRIVIAL-GARBAGE {124BF099}> as TRIVIAL-GARBAGE ; registering #<SYSTEM TRIVIAL-GARBAGE-TESTS {1272DDF1}> as ; TRIVIAL-GARBAGE-TESTS ; loading system definition from ; /Users/malc/opt/lisp/sbcl/sbcl-inst-1.0.23-x86-darwin/lib/sbcl/site-systems/cffi.asd ; into #<PACKAGE "ASDF0"> ; registering #<SYSTEM CFFI {119164D1}> as CFFI ; loading system definition from ; /Users/malc/opt/lisp/sbcl/sbcl-inst-1.0.23-x86-darwin/lib/sbcl/site-systems/babel.asd ; into #<PACKAGE "ASDF0"> ; registering #<SYSTEM BABEL {11CA1489}> as BABEL ; loading system definition from ; /Users/malc/opt/lisp/sbcl/sbcl-inst-1.0.23-x86-darwin/lib/sbcl/site-systems/alexandria.asd ; into #<PACKAGE "ASDF0"> ; registering #<SYSTEM :ALEXANDRIA {1276FF09}> as ALEXANDRIA ; loading system definition from ; /Users/malc/opt/lisp/sbcl/sbcl-inst-1.0.23-x86-darwin/lib/sbcl/site-systems/trivial-features.asd ; into #<PACKAGE "ASDF0"> ; registering #<SYSTEM TRIVIAL-FEATURES {129ED569}> as TRIVIAL-FEATURES NIL CL-USER> (asdf:operate 'asdf:load-op :gsll-tests) ; loading system definition from ; /Users/malc/opt/lisp/asdf-registry/gsll-tests.asd into #<PACKAGE "ASDF0"> ; registering #<SYSTEM "gsll-tests" {11B966E9}> as gsll-tests ; loading system definition from ; /Users/malc/opt/lisp/asdf-registry/lisp-unit.asd into #<PACKAGE "ASDF0"> ; registering #<SYSTEM :LISP-UNIT {12484CB9}> as LISP-UNIT NIL CL-USER> (in-package :gsll) #<PACKAGE "GSLL"> GSL> (lisp-unit:assert-numerical-equal '(1 2 3) '(1 2 4)) '(1 2 4) failed: Expected (1 2 3) but saw (1 2 4) NIL GSL> (lisp-unit:assert-numerical-equal (make-marray 'double-float :initial-contents '(1 2 3)) (make-marray 'double-float :initial-contents '(1 2 4))) T
Can you guys reproduce this? I'm SBCL 1.0.23 on OS X if that makes a difference..
Cheers
Malcolm
-
Liam Healy -
Malcolm Reynolds