In today’s Programming Praxis exercise we have to do some numerical integration. Let’s get started.

Since the first three algorithms all work roughly the same, and as a programmer I hate repeating myself, we’re going to make a generic integration function to abstract out the common stuff.

int combine f a b n = sum $ map g [0..n - 1] where
w = (b - a) / n
lo i = a + w * i
g i = w * combine (f $ lo i) (f $ lo i + w/2) (f $ lo i + w)

Now we can just define the three simple integration methods in terms of what we want to do with the low, mid and high points.

intRect = int (\_ m _ -> m)
intTrap = int (\l _ h -> (l + h) / 2)
intSimp = int (\l m h -> (l + 4 * m + h) / 6)

The intAdapt function is pretty much the same as the Scheme version.

intAdapt m f a b epsilon = if abs (g10 - g5) < e then g10 else
intAdapt m f a mid (Just e) + intAdapt m f mid b (Just e)
where g5 = m f a b 5
g10 = m f a b 10
mid = (a + b) / 2
e = maybe 1e-7 id epsilon

Using adaptive integration for prime counting:

approxPi n = round $ intAdapt intSimp (recip . log) 2 n Nothing

And finally a test to show that everything’s working properly.

main = do print $ intRect cube 0 1 10000
print $ intTrap cube 0 1 10000
print $ intSimp cube 0 1 10000
print $ approxPi 1e21
where cube x = x * x * x

### Like this:

Like Loading...

*Related*

Tags: bonsai, code, Haskell, integration, kata, numerical, praxis, programming

This entry was posted on February 9, 2010 at 10:28 am and is filed under Programming Praxis. You can follow any responses to this entry through the RSS 2.0 feed.
You can leave a response, or trackback from your own site.

## Leave a Reply