Rosetta Code

Numerical integration: Difference between revisions

Page Discussion

Numerical integration (view source)

Revision as of 13:22, 6 March 2017

282 bytes removed ,  7 years ago
→‎{{header|Haskell}}: Added section labels to output, pruned some redundancy in output code
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 f xs ws = sum [w * f x | (x,w) <- zip xs ws]
:: Fractional a
=> (a1 -> a) -> [a1] -> [a] -> a
approx f xs ws =
sum
[ w * f x
| (x, w) <- zip xs ws ]
 
integrateOpen
integrateOpen :: Fractional a => a -> [a] -> (a -> a) -> a -> a -> Int -> a
:: Fractional a
integrateOpen v vs f a b n = approx f xs ws * h / v where
=> a -> [a] -> (a -> a) -> a -> a -> Int -> a
m = fromIntegral (length vs) * n
integrateOpen v vs f a b n = approx f xs ws * h / v
h = (b-a) / fromIntegral m
where
ws = concat $ replicate n vs
m = fromIntegral (length vs) * n
c = a + h/2
xs = [c + h *= fromIntegral(b i- |a) i/ <-fromIntegral [0..m-1]]
ws = concat $ replicate n vs
c = a + h / 2
xs =
[ c + h * fromIntegral i
| i <- [0 .. m - 1] ]
 
integrateClosed
integrateClosed :: Fractional a => a -> [a] -> (a -> a) -> a -> a -> Int -> a
:: Fractional a
integrateClosed v vs f a b n = approx f xs ws * h / v where
=> a -> [a] -> (a -> a) -> a -> a -> Int -> a
m = fromIntegral (length vs - 1) * n
integrateClosed v vs f a b n = approx f xs ws * h / v
h = (b-a) / fromIntegral m
where
ws = overlap n vs
m = fromIntegral (length vs - 1) * n
xs = [a + h * fromIntegral i | i <- [0..m]]
h = (b - a) / fromIntegral m
ws = overlap n vs
xs =
[ a + h * fromIntegral i
| i <- [0 .. m] ]
 
overlap
overlap :: Num a => Int -> [a] -> [a]
:: Num a
overlap n [] = []
=> Int -> [a] -> [a]
overlap n (x:xs) = x : inter n xs where
overlap n [] = []
inter 1 ys = ys
interoverlap n [] (x:xs) = x : inter (n-1) xs
where
inter n [y] = (x+y) : inter (n-1) xs
inter n1 (y:ys) = y : inter n ys
inter n [] = x : inter (n - 1) xs
 
inter n [y] = (x + y) : inter (n - 1) xs
intLeftRect = integrateClosed 1 [1,0]
inter n (y:ys) = y : inter n ys
intMidRect = integrateOpen 1 [1]
intRightRect = integrateClosed 1 [0,1]
intTrapezium = integrateClosed 2 [1,1]
intSimpson = integrateClosed 3 [1,4,1]
 
uncurry4 :: (t1 -> t2 -> t3 -> t4 -> t) -> (t1, t2, t3, t4) -> t
uncurry4 f ~(a, b, c, d) = f a b c d
 
-- TEST ----------------------------------------------------------------------
main = do
ms
let m1 = "rectangular left: "
:: Fractional a
let m2 = "rectangular middle: "
=> [(String, (a -> a) -> a -> a -> Int -> a)]
let m3 = "rectangular right: "
ms =
let m4 = "trapezium: "
[ ("rectangular left", integrateClosed 1 [1, 0])
let m5 = "simpson: "
, ("rectangular middle", integrateOpen 1 [1])
 
, ("rectangular right", integrateClosed 1 [0, 1])
let arg1 = ((\x -> x ^ 3), 0, 1, 100)
, ("trapezium", integrateClosed 2 [1, 1])
putStrLn $ m1 ++ (show $ uncurry4 intLeftRect arg1)
, ("simpson", integrateClosed 3 [1, 4, 1])
putStrLn $ m2 ++ (show $ uncurry4 intMidRect arg1)
]
putStrLn $ m3 ++ (show $ uncurry4 intRightRect arg1)
putStrLn $ m4 ++ (show $ uncurry4 intTrapezium arg1)
putStrLn $ m5 ++ (show $ uncurry4 intSimpson arg1)
putStrLn ""
 
let arg2 = ((\x -> 1 / x), 1, 100, 1000)
putStrLn $ m1 ++ (show $ uncurry4 intLeftRect arg2)
putStrLn $ m2 ++ (show $ uncurry4 intMidRect arg2)
putStrLn $ m3 ++ (show $ uncurry4 intRightRect arg2)
putStrLn $ m4 ++ (show $ uncurry4 intTrapezium arg2)
putStrLn $ m5 ++ (show $ uncurry4 intSimpson arg2)
putStrLn ""
 
let arg3 = ((\x -> x), 0, 5000, 5000000)
putStrLn $ m1 ++ (show $ uncurry4 intLeftRect arg3)
putStrLn $ m2 ++ (show $ uncurry4 intMidRect arg3)
putStrLn $ m3 ++ (show $ uncurry4 intRightRect arg3)
putStrLn $ m4 ++ (show $ uncurry4 intTrapezium arg3)
putStrLn $ m5 ++ (show $ uncurry4 intSimpson arg3)
putStrLn ""
 
expressions
let arg4 = ((\x -> x), 0, 6000, 6000000)
:: (Fractional a, Num t, Num t1, Num t2)
putStrLn $ m1 ++ (show $ uncurry4 intLeftRect arg4)
=> [(String, (a -> a, t, t1, t2))]
putStrLn $ m2 ++ (show $ uncurry4 intMidRect arg4)
expressions =
putStrLn $ m3 ++ (show $ uncurry4 intRightRect arg4)
[ ("x^3", ((^ 3), 0, 1, 100))
putStrLn $ m4 ++ (show $ uncurry4 intTrapezium arg4)
, ("1/x", ((1 /), 1, 100, 1000))
putStrLn $ m5 ++ (show $ uncurry4 intSimpson arg4)</lang>
, ("x", (id, 0, 5000, 500000))
, ("x", (id, 0, 6000, 600000))
]
 
main :: IO ()
Output:
main =
<pre>rectangular left: 0.24502500000000005
mapM_
rectangular middle: 0.24998750000000006
(\(s, e@(f, a, b, n)) -> do
rectangular right: 0.25502500000000006
putStrLn $
trapezium: 0.25002500000000005
simpson: justifyLeft 20 ' ' 0.25000000000000006("f(x) = " ++ s) ++
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 left: 4.65499105751468
rectangular middleleft: 4.60476254867837665499105751468
rectangular rightmiddle: 4.55698105751468604762548678376
rectangular right: 4.55698105751468
trapezium: 4.605986057514681
simpsontrapezium: 4.605170384957134605986057514681
simpson: 4.605170384957135
 
f(x) = x 0.0 5000.0 500000
rectangular left: 1.24999975e7
rectangular middleleft: 1.25e72499975000000006e7
rectangular rightmiddle: 1.25000025e72499999999999993e7
rectangular right: 1.2500025000000006e7
trapezium: 1.25e7
simpsontrapezium: 1.2499999999999993e72500000000000006e7
simpson: 1.2499999999999998e7
 
f(x) = x 0.0 6000.0 600000
rectangular left: 1.7999997000000004e7
rectangular middleleft: 1.7999999999999993e77999970000000004e7
rectangular rightmiddle: 1.8000003000000004e77999999999999993e7
rectangular right: 1.8000030000000004e7
trapezium: 1.8000000000000004e7
simpsontrapezium: 1.7999999999999993e7</pre>8000000000000004e7
simpson: 1.7999999999999996e7</pre>
Runtime: about 7 seconds.
 
Hout
9,659

edits

Retrieved from "https://rosettacode.org/wiki/Special:MobileDiff/52526"