[Haskell-cafe] Gauss Elimination -> More Clean2Haskell
Philippos Apolinarius
1 --- The program must be implemented using arrays.� Update must be done in place, in order to minimize use of the garbage collector. I have used Data.Array.IO, but I guess that� Data.Array.ST is better. Is it easy to rewrite the program in order to use DataArray.ST?
2 --� I liked very much the forM_ monad. I wonder if there is an accumulating monad that play the role of a program like leftSide.
3 -- Clean program almost double its speed, if one uses strictness annotations. Is it possible to use similar anotations for Haskell?
Here is how I compiled the program:
ghc -O2 gel.hs --make
In order to run the program with 2000 equations,� type
gel.exe 2000 +RTS -sstderr
The program will create a linear system with 2000 equations so that all elements of the solution is equal to 1. It prints 20 elements of the solution.
Here is the program:
{- File: gel.hs
Compilation: ghc -O2 gel.hs --make
Run: time gel.exe 2000
-}
import Control.Monad
import Data.Array.IO
import System.IO
import System.Random
import System (getArgs)
prtSol i n1 arr | i < 1= return ()
prtSol i n1 arr= do
��� b <- readArray arr (i, n1)
��� putStr ((show b)++" ")
��� prtSol (i-1) n1 arr
fillArray xs s (i, n) (j, m) arr |� i > n= return ()
fillArray xs s (i,n) (j, m) arr | i==n && j>m= do
� writeArray arr (i, j) s
� return ()
fillArray xs s (i, n) (j, m) arr | j > m� = do
�� writeArray arr (i, j) s
�� fillArray xs 0.0 (i+1, n) (1, m) arr
fillArray (val:xs) s (i, n) (j, m) arr= do
�� writeArray arr (i, j) val
�� fillArray xs (s+val) (i, n) (j+1, m) arr
sLU arr n= sIJ 2 1 2 n arr
sIJ i j k n arr | i > n = return ()
sIJ i j k n arr | k > n = sIJ (i+1) i (i+1) n arr
sIJ i j k n arr = do
� im <- pmax (j+1) j
� swap j im 1
� a <- readArray arr (k, j)
� forM_ [j..n+1] $� \l -> do
����� ajj <- readArray arr (j, j)
����� ajl <- readArray arr (j, l)
����� akl <- readArray arr (k, l)
����� writeArray arr (k, l) (akl-a*(ajl/ajj))
� sIJ i j (k+1) n arr where
���� pmax line imax | line > n = return imax
���� pmax line imax = do
������ alj <- readArray arr (line, j)
������ aij <- readArray arr (imax, j)
������ if (abs alj)> (abs aij)
��������� then pmax (line+1) line
��������� else pmax (line+1) imax
���� swap r s q | q>n+1 = return ()
���� swap r s q | r==s = return ()
���� swap r s q = do
������� arq <- readArray arr (r,q)
������� asq <- readArray arr (s,q)
������� writeArray arr (s,q) arq
������� writeArray arr (r,q) asq
������� swap r s (q+1)
����
�
leftSide acc i j n arr | j>n= return acc
leftSide acc i j n arr = do
�� v <- readArray arr (j, n+1)
�� a <- readArray arr (i, j)
�� leftSide (acc-v*a) i (j+1) n arr
solvit i n arr | i<1= return ()
solvit i n arr= do
�� a <- readArray arr (i, i)
�� acc <- readArray arr (i, n+1)
�� v <- leftSide acc i (i+1) n arr
�� writeArray arr (i, n+1) (v/a)
�� solvit (i-1) n arr
rnList :: (Double, Double) -> IO [Double]
rnList r=getStdGen>>=(\x->return(randomRs r x))
dims [input] = (read input, read input)
dims _ = (1000, 1000)
main = do
���� xs <- rnList (1.0,1000.0)
���� args <- getArgs
���� let (n, m)= dims args
���� arr <- newArray_ ((1,1),(n,m+1)) ::
��������������� IO (IOUArray (Int, Int) Double)
���� fillArray xs 0.0 (1,n) (1,m) arr
���� sLU arr n
���� solvit n n arr
���� prtSol (min 20 n) (n+1) arr
���� print "Done"
__________________________________________________________________
Be smarter than spam. See how smart SpamGuard is at giving junk email the boot with the All-new Yahoo! Mail. Click on Options in Mail and switch to New Mail today or register for free at http://mail.yahoo.ca
Ryan Ingram
<phi5...@yahoo.ca> wrote:
> 1 --- The program must be implemented using arrays.
> Update must be done in place, in order to minimize use of the garbage collector.
> I have used Data.Array.IO, but I guess that� Data.Array.ST is better. Is it easy to
> rewrite the program in order to use DataArray.ST?
It should be pretty easy as long as the rest is pure; ST is a subset
of I/O that deals with algorithms that have mutable variables/arrays
but no observable side-effects. ST safely guarantees that no
side-effects escape to somewhere they can be observed through a clever
type-system trick.
> 2 --� I liked very much the forM_ monad. I wonder if there is an accumulating monad that play the role of a program like leftSide.
forM_ is just a function; it works for all monads. Probably just a
terminology error?
> 3 -- Clean program almost double its speed, if one uses strictness annotations. Is it possible to use similar anotations for Haskell?
Yes. The common ways to do this are to use ! annotations in data
structures, like so:
] data Foo s = Foo !Int !Int !(STArray s (Int,Int) Double)
You also can use seq and/or $! to guide the evaluation order of your
expressions:
x <- readArray a (1,1)
writeArray a (1,1) $! (x+1) -- forces x+1 to evaluate instead of
writing a thunk.
If you really care about speed, you probably want to look into
STUArrays; these store unboxed values and should be about as fast as a
C array.
Now to the stylistic comments:
You can use guards better to not repeat yourself so often:
> prtSol i n1 arr | i < 1= return ()
> prtSol i n1 arr= do
> ��� b <- readArray arr (i, n1)
> ��� putStr ((show b)++" ")
> ��� prtSol (i-1) n1 arr
becomes
] prtSol i n1 arr
] | i < 1 = return ()
] | otherwise = do
] b <- readArray arr (i, n1)
] putStr ((show b)++" ")
] prtSol (i-1) n1 arr
Similarily:
> fillArray xs s (i, n) (j, m) arr |� i > n= return ()
> fillArray xs s (i,n) (j, m) arr | i==n && j>m= do
this branch doesn't need "do" because writeArray returns ()
> � writeArray arr (i, j) s
> � return ()
> fillArray xs s (i, n) (j, m) arr | j > m� = do
> �� writeArray arr (i, j) s
> �� fillArray xs 0.0 (i+1, n) (1, m) arr
> fillArray (val:xs) s (i, n) (j, m) arr= do
> �� writeArray arr (i, j) val
> �� fillArray xs (s+val) (i, n) (j+1, m) arr
] fillArray xs s (i, n) (j, m) arr
] | i > n= return ()
] | i==n && j>m = writeArray arr (i, j) s
] | j > m = do
] writeArray arr (i, j) s
] fillArray xs 0.0 (i+1, n) (1, m) arr
] | otherwise = do
] writeArray arr (i, j) val
] fillArray xs (s+val) (i, n) (j+1, m) arr
I'll let someone else show you how to build this into a fold.
-- ryan
_______________________________________________
Haskell-Cafe mailing list
Haskel...@haskell.org
http://www.haskell.org/mailman/listinfo/haskell-cafe
Don Stewart
> import Control.Monad
> import Data.Array.IO
You're using boxed arrays here -- try Data.Array.Unboxed and IOUArray's
Felipe Lessa
> ] prtSol i n1 arr
> ] | i < 1 = return ()
> ] | otherwise = do
> ] b <- readArray arr (i, n1)
> ] putStr ((show b)++" ")
> ] prtSol (i-1) n1 arr
Which is just
> prtSol i n1 arr = unless (i < 1) $ do
> b <- readArray arr (i, n1)
> putStr ((show b)++" ")
> prtSol (i-1) n1 arr
--
Felipe.