
Hello Martin, Sunday, March 12, 2006, 8:49:15 PM, you wrote: MP> 1. Haskell-nativeness: I have had some issues compiling and linking with gsl MP> libraries on 64-bit platforms. Also, it would be quite interesting to gauge MP> haskell's effectiveness as a scientific computing platform, in particular MP> making use unboxed arrays as representation. it's a bad idea. simple test of vector addition shows the 20x performance loss comparing GHC and GCC. so you can consider this as learning example but not as a real-world lib. btw, i've proposed making changes in ghc that will change this situation, at least for simple loops. but that is only words for this moment. if you want to allow ghc became the really good choice for writing matrix libraries, it's better to take part in ghc's own development :) MP> 2. Use of monads: in the afore-mentioned libraries, monads are either ignored, MP> or are restricted to the IO monad. IO monad is just ST monad specialized to RealWorld state. you can easily use this libs in ST monad with help of unsafeIOToST. surprised? :) you will laugh even more if i say that STUArray in Hugs implemented in just this way - by using peek/poke IO operations wrapped in unsafeIOToST. just small excerpt from this implementation: specialIOToST :: IO a -> ST s a specialIOToST = unsafeCoerce type BytePtr = ForeignPtr Word8 data MutableByteArray s = MutableByteArray !Int !BytePtr newMutableByteArray :: Int -> ST s (MutableByteArray s) newMutableByteArray size = do fp <- specialIOToST (mallocForeignPtrBytes size) return (MutableByteArray size fp) readMutableByteArray :: Storable e => MutableByteArray s -> Int -> ST s e readMutableByteArray (MutableByteArray _ fp) i = specialIOToST $ withForeignPtr fp $ \a -> peekElemOff (castPtr a) i MP> I am taking a different approach, starting MP> my code from the ST monad, and eventually this will be generalized to work with MP> the IO monad. I think the ST monad is a good monad to be able to perform MP> computations on matrices that update them, allowing efficient, in-place MP> algorithms, but it also provides the benifit of not being a one-way street like MP> the IO monad is. Being a relative newcomer to haskell, I would be interested to MP> hear any thoughts as to whether this is a good/bad idea. it's a right way (and even obvious way to one who knows how this all work). btw, i developed general i/o and serialization library that works both in IO and ST monad, and now i almost finished rewriting of arrays and references library, which internally uses the same concept of providing general code that works in both monads. just a small excerpt from my code: -- | Unboxed mutable arrays data UnboxedMutableArray s i e = UMA !i !i !(MUVec s e) instance (STorIO m s) => HasMutableBounds (UnboxedMutableArray s) m where getBounds (UMA l u _) = return (l,u) instance (STorIO m s, Unboxed e) => MArray (UnboxedMutableArray s) e m where newArray_ (l,u) = do arr <- allocUnboxed (rangeSize (l,u)) return (UMA l u arr) unsafeRead (UMA _ _ arr) index = readUnboxed arr index unsafeWrite (UMA _ _ arr) index = writeUnboxed arr index -- | Unboxed mutable arrays in ST monad type STUArray = UnboxedMutableArray -- | Unboxed mutable arrays in IO monad type IOUArray = IOSpecific3 UnboxedMutableArray my library also provides monad-independent references. i.e. you can write monadic code that works with references and this code can be runned without any problems both in IO and ST monads: -- This section demonstrates running of monad-independent algorithm -- `test_Ref` in IO and ST monads test_Ref 3 >>= print print $ runST (test_Ref 4) i can also add support for monad-independent array manipulations so that you can write code that will work in both monads. it was my old idea but i had no clients for it :) MP> Now to my question: I would like to represent a matrix as a wrapper around a MP> block, which in turn is just an unboxed array. Here are the definitions for a MP> matrix in ST and outside of a monad, respectively: MP> type MBlock s = STUArray s Int Double MP> data MMatrix s = MMatrix Int Int (MBlock s) MP> type Block = UArray Int Double MP> data Matrix = Matrix Int Int Block seems that you don't know that Haskell's indexes can be a tuples :) type Matrix a = UArray (Int,Int) a type Matrix3d a = UArray (Int,Int,Int) a MP> Now, I have started by providing some fairly low-level routines in the ST MP> monad. My problem is depressingly simple: I would like to retrieve a matrix MP> from the ST monad for use outside of it (in particular, to pretty-print it). MP> Now, this is easy for STUArrays (just use runSTUArray), but I'm not sure of how MP> to do it for a type that *encloses* an STUArray. may be you need to look inside runSTUArray's definition? ;) that is from one doc i wrote: "Operations that creates/updates immutable arrays just creates them as mutable arrays in ST monad, make all required updates on this array and then use unsafeFreeze before returning array from runST." so, you can return anything else together with array returned by unsafeFreeze. of course, because you should use tuples for indexing, it has only theoretical interest. but i strongly recommend you to read entire data.array.* sources to know about all intrinsics of arrays library and in particular freeze/thaw tale MP> runSTMatrix :: ST s (Matrix s) -> Matrix runSTMatrix :: ST s (MMatrix s) -> Matrix runSTMatrix a = runST ( do (MMatrix i j mblock) <- a block <- unsafeFreeze mblock return (Matrix i j block) ) -- Best regards, Bulat mailto:Bulat.Ziganshin@gmail.com