Re: haskell blas bindings: does iomatrix gemv transposing of matrix a?
Hi Anatoly, I made the decision to make "herm" an O(1) operation. This means you don't have to pass transpose arguments to the multiplication routines. When you do, for example:
let a = listMatrix (2,3) [1..6] x = listVector 2 [1, -1] in herm a <*> x
this gets implemented as a call to gemv with transa set to "ConjTrans".
BLAS actually supports "NoTrans", "Trans", and "ConjTrans". For real data types, "Trans" is the same as "ConjTrans". For complex data types, they are different. Unfortunately, it is impossible to call a routine with "Trans" in the blas bindings. I explicitly did not want people to have to pass extra parameters to indicate transpose/conjugacy information, and I wanted "conjugate transpose" to be an O(1) operation. Ideally, "transpose" would be O(1), as well. However, since BLAS doesn't have an option to conjugate the input matrix but not transpose it, only one of "herm" and "trans" can be O(1).
I think in hmatrix, "trans" is O(1), but not "herm".
Patrick
p.s. I've been busy with other work lately, and I probably won't make a release soon, but I did a pretty big overhaul of the library to support operations in the ST monad. For impure operations, the functions of all the operations are mostly the same, but the types are all different. (For pure operations, there are only minor changes.) If you want the development version, you can get the darcs repo at http://stat.stanford.edu/~patperry/code/blas .
The new version is pretty typeclass-heavy, since that's the only way I know how to support both ST and IO. Consequently, there have been some performance regressions. I have some optimization ideas (in the "TODO") file, but I do not have time to implement them right now. If you or anyone else would like to help with this or anything else, I would be glad to have you aboard the development "team".
----- Original Message -----
From: "Anatoly Yakovenko"
I made the decision to make "herm" an O(1) operation. This means you don't have to pass transpose arguments to the multiplication routines. When you do, for example:
let a = listMatrix (2,3) [1..6] x = listVector 2 [1, -1] in herm a <*> x
this gets implemented as a call to gemv with transa set to "ConjTrans".
Ah, i see, i didn't see that function. That's pretty slick actually.
The new version is pretty typeclass-heavy, since that's the only way I know how to support both ST and IO. Consequently, there have been some performance regressions. I have some optimization ideas (in the "TODO") file, but I do not have time to implement them right now. If you or anyone else would like to help with this or anything else, I would be glad to have you aboard the development "team".
I would be glad to help. Its probably about time i learned the type system anyways. Anatoly
is there anyway the modifyWith functions could work on uboxed types?
aeyakovenko:
is there anyway the modifyWith functions could work on uboxed types?
If they're inlined, the modify functions on boxed types may well end up unboxed. What's the particular problem you're having? -- Don
is there anyway the modifyWith functions could work on uboxed types?
If they're inlined, the modify functions on boxed types may well end up unboxed.
What's the particular problem you're having?
well, after inspecting a little further its not so bad actually. i was comparing module Main where import qualified Data.Vector.Dense.IO as Vector import Control.Monad e = exp 1.0 sigmoid xx = 1.0 / (1 + (e ** (1.0 * xx))) type Vec = Vector.IOVector Int Double main = do let size = 100000 ff::Vec <- Vector.newListVector size $ repeat 0.5 replicateM_ 1000 $ Vector.modifyWith (sigmoid) ff putStrLn $ "done" to this: #include "math.h" #include "stdlib.h" #include "stdio.h" double sigmoid(double xx) { return 1.0 / (1.0 + (pow(M_E, (1.0 * xx)))); } int main() { int size = 100000; int times = 1000; int ii,jj; double* array = malloc(sizeof(double)*size); for(jj = 0; jj < size; ++jj) { array[jj] = 0.5; } for(ii = 0; ii < times; ++ii) { for(jj = 0; jj < size; ++jj) { array[jj] = sigmoid(array[jj]); } } printf("done\n"); } the haskell version does ok, or 0m37.937s vs 0m23.492s in C. I am compiling with these options: -O2 -fexcess-precision -funbox-strict-fields -fglasgow-exts -fbang-patterns -prof -auto-all, and O2 for gcc.
On Wed, 24 Sep 2008, Anatoly Yakovenko wrote:
is there anyway the modifyWith functions could work on uboxed types?
If they're inlined, the modify functions on boxed types may well end up unboxed.
What's the particular problem you're having?
well, after inspecting a little further its not so bad actually. i was comparing
module Main where
import qualified Data.Vector.Dense.IO as Vector import Control.Monad
e = exp 1.0 sigmoid xx = 1.0 / (1 + (e ** (1.0 * xx)))
I think 'exp' is generally more efficient and less ambiguous than 'e **'.
Anatoly Yakovenko wrote:
e = exp 1.0 sigmoid xx = 1.0 / (1 + (e ** (1.0 * xx)))
That 1.0 * xx caught my eye. In case this was an oversight on your part: if you mean the usual sigmoid function, that should be 1.0 / (1 + (e ** (0.0 - x))). -- Joe Buehler
e = exp 1.0 sigmoid xx = 1.0 / (1 + (e ** (1.0 * xx)))
That 1.0 * xx caught my eye.
In case this was an oversight on your part: if you mean the usual sigmoid function, that should be 1.0 / (1 + (e ** (0.0 - x))).
i had a different constant there before.
participants (5)
-
Anatoly Yakovenko -
Don Stewart -
Henning Thielemann -
Joe Buehler -
Patrick O'Regan Perry