Numerical integration: Difference between revisions
Content deleted Content added
SqrtNegInf (talk | contribs) m →{{header|Perl 6}}: cosmetic fix |
→{{header|Haskell}}: Added section labels to output, pruned some redundancy in output code |
||
Line 2,276:
The whole program:
<lang haskell>approx
:: Fractional a
=> (a1 -> a) -> [a1] -> [a] -> a
approx f xs ws =
sum
[ w * f x
| (x, w) <- zip xs ws ]
integrateOpen
:: Fractional a
=> a -> [a] -> (a -> a) -> a -> a -> Int -> a
integrateOpen v vs f a b n = approx f xs ws * h / v
where
m = fromIntegral (length vs) * n
ws = concat $ replicate n vs
c = a + h / 2
xs =
[ c + h * fromIntegral i
| i <- [0 .. m - 1] ]
integrateClosed
:: Fractional a
=> a -> [a] -> (a -> a) -> a -> a -> Int -> a
integrateClosed v vs f a b n = approx f xs ws * h / v
where
m = fromIntegral (length vs - 1) * n
h = (b - a) / fromIntegral m
ws = overlap n vs
xs =
[ a + h * fromIntegral i
| i <- [0 .. m] ]
overlap
:: Num a
=> Int -> [a] -> [a]
overlap n [] = []
where
inter
inter n [] = x : inter (n - 1) xs
inter n [y] = (x + y) : inter (n - 1) xs
inter n (y:ys) = y : inter n ys
uncurry4 :: (t1 -> t2 -> t3 -> t4 -> t) -> (t1, t2, t3, t4) -> t
uncurry4 f ~(a, b, c, d) = f a b c d
-- TEST ----------------------------------------------------------------------
ms
:: Fractional a
=> [(String, (a -> a) -> a -> a -> Int -> a)]
ms =
[ ("rectangular left", integrateClosed 1 [1, 0])
, ("rectangular middle", integrateOpen 1 [1])
, ("rectangular right", integrateClosed 1 [0, 1])
, ("trapezium", integrateClosed 2 [1, 1])
, ("simpson", integrateClosed 3 [1, 4, 1])
]
expressions
:: (Fractional a, Num t, Num t1, Num t2)
=> [(String, (a -> a, t, t1, t2))]
expressions =
[ ("x^3", ((^ 3), 0, 1, 100))
, ("1/x", ((1 /), 1, 100, 1000))
, ("x", (id, 0, 5000, 500000))
, ("x", (id, 0, 6000, 600000))
]
main :: IO ()
main =
mapM_
(\(s, e@(f, a, b, n)) -> do
putStrLn $
unwords (show <$> [a, b]) ++ " " ++ show n
mapM_
(\(s, integration) ->
putStrLn
(justifyLeft 20 ' ' (s ++ ":") ++ show (uncurry4 integration e)))
ms
putStrLn [])
expressions
where
justifyLeft n c s = take n (s ++ replicate n c)</lang>
{{Out}}
<pre>f(x) = x^3 0.0 1.0 100
rectangular left: 0.24502500000000005
rectangular middle: 0.24998750000000006
rectangular right: 0.25502500000000006
trapezium: 0.25002500000000005
simpson: 0.25000000000000006
f(x) = 1/x 1.0 100.0 1000
rectangular
rectangular
rectangular right: 4.55698105751468
simpson: 4.605170384957135
f(x) = x 0.0 5000.0 500000
rectangular
rectangular
rectangular right: 1.2500025000000006e7
simpson: 1.2499999999999998e7
f(x) = x 0.0 6000.0 600000
rectangular
rectangular
rectangular right: 1.8000030000000004e7
simpson: 1.7999999999999996e7</pre>
Runtime: about 7 seconds.
|