I've been dabbling with implementing the "n-body" benchmark over at
the "Great Computer Language Shootout"...http://shootout.alioth.debian.org/benchmark.php?test=nbody ...and I've stumbled on to a problem with space leaks. The code below works great for low numbers of iterations, but crashes and burns with large numbers. You can temporarily fix it up by specifying larger stack and heap sizes (with the +RTS -K & -H options), but this only gets you so far. I haven't been able to glean much information from the profiling reports. You can see the reports I generated (like the biography report) over at <http://sleepingsquirrel.org/nbody>. I've tried several different variations on the "advance" function, as well as trying "seq" in numerous spots, but I haven't stumbled on to the right combination yet. Any pointers? Thanks, Greg Buchholz >-- n-body simulation >-- compile: ghc -O2 -o nbody nbody.hs >-- run: nbody 100000 > >import System(getArgs) >import Numeric > >data Vec = V !Double !Double !Double deriving Show >instance Num Vec where > (V a b c) + (V x y z) = (V (a+x) (b+y) (c+z)) > (V a b c) - (V x y z) = (V (a-x) (b-y) (c-z)) > fromInteger 0 = V 0.0 0.0 0.0 -- for sum function >instance Eq Vec where (V a b c) == (V x y z) = (a==x) && (b==y) && (c==z) > >dot (V a b c) (V x y z) = a*x + b*y + c*z >scale (V a b c) n = V (n*a) (n*b) (n*c) >mag (V x y z) = sqrt (x*x + y*y + z*z) > >data Planet = Planet Vec Vec Double deriving Show --Position Velocity Mass >dist (Planet p1 _ _) (Planet p2 _ _) = mag $ p1 - p2 >mass (Planet _ _ m) = m >vel (Planet _ v _) = v >pos (Planet p _ _) = p > >main = do > [arg] <- getArgs > let iter = read arg > let n = 5::Int > let bodies = offset_momentum n [sun, jupiter, saturn, neptune, uranus] > let begin = energy n bodies > putStrLn $ showFFloat (Just 9) begin "" > let final = (iterate (advance n 0.01) bodies) > let end = energy n (final !! iter) > putStrLn $ showFFloat (Just 9) end "" > >days_per_year = 365.24 >solar_mass = 4.0 * pi * pi > >advance :: Int -> Double -> [Planet] -> [Planet] >advance n dt (p:ps) = take n ((Planet (pos p + delta_x) > (new_v) (mass p) ) : advance n dt (ps++[p])) > where > delta_v = sum (map (\q -> > (pos p - pos q) `scale` ((mass q)*dt/(dist p q)^3)) ps) > new_v = (vel p) - delta_v > delta_x = new_v `scale` dt > >energy:: Int -> [Planet] -> Double >energy n ps = kinetic - potential > where > kinetic = 0.5 * (sum (map (\q->(mass q)*((vel q) `dot` (vel q))) ps)) > potential = sum [(mass (ps!!i))*(mass (ps!!j))/(dist (ps!!i) (ps!!j)) > | i<-[0..n-1], j<-[i+1..n-1]] > >offset_momentum n ((Planet p v m):ps) = (Planet p new_v m):ps > where new_v = (sum (map (\n->(vel n) `scale` (mass n)) ps)) > `scale` ((-1.0)/solar_mass) > >{- >rotations 0 _ = [] >rotations n xss@(x:xs) = flipped:rotations (n-1) flipped > where flipped = xs++[x] > >advance :: Int -> Double -> [Planet] -> [Planet] >advance n dt pss = map (adv dt) (rotations n pss) > >adv dt (p:ps) = Planet (pos p + delta_x) (new_v) (mass p) > where > delta_v = sum (map (\q -> > (pos p - pos q) `scale` ((mass q)*dt/(dist p q)^3)) ps) > new_v = (vel p) - delta_v > delta_x = new_v `scale` dt >--} >jupiter = (Planet > (V 4.84143144246472090e+00 (-1.16032004402742839e+00) > (-1.03622044471123109e-01)) > (V ( 1.66007664274403694e-03 * days_per_year) > ( 7.69901118419740425e-03 * days_per_year) > ((-6.90460016972063023e-05) * days_per_year)) > (9.54791938424326609e-04 * solar_mass)) > >saturn = (Planet > (V 8.34336671824457987e+00 4.12479856412430479e+00 > (-4.03523417114321381e-01)) > (V (-2.76742510726862411e-03 * days_per_year) > (4.99852801234917238e-03 * days_per_year) > (2.30417297573763929e-05 * days_per_year)) > (2.85885980666130812e-04 * solar_mass)) > >uranus = (Planet > (V 1.28943695621391310e+01 (-1.51111514016986312e+01) > (-2.23307578892655734e-01)) > (V (2.96460137564761618e-03 * days_per_year) > (2.37847173959480950e-03 * days_per_year) > (-2.96589568540237556e-05 * days_per_year)) > (4.36624404335156298e-05 * solar_mass)) > >neptune = (Planet > (V 1.53796971148509165e+01 (-2.59193146099879641e+01) > 1.79258772950371181e-01) > (V (2.68067772490389322e-03 * days_per_year) > (1.62824170038242295e-03 * days_per_year) > (-9.51592254519715870e-05 * days_per_year)) > (5.15138902046611451e-05 * solar_mass)) > >sun = (Planet (V 0.0 0.0 0.0) (V 0.0 0.0 0.0) solar_mass) _______________________________________________ Haskell-Cafe mailing list [email protected] http://www.haskell.org/mailman/listinfo/haskell-cafe
