The most compact algorithm for Simpson's rule??

[Dieter]

(12-12-2015 04:24 PM)Namir Wrote:  Here is a pseudo code for what I hope to be the most compact way to implement Simpson's rule:

I'm not sure if this counts as "more compact", but the alternating factors k=4, 2, 4, 2,... can be generated by subtraction from the sum of both, i.e. k := 6 – k.

I would also use a simple loop counter instead of continuously adding h (which may lead to roundoff errors).

And finally, the initial value of "sum" of course has to be set to the sum of f(a) and f(b), not the difference.

Edit: I now see this is not an error since you calculate f(b) twice, once before the loop and once within it.
The latter value is added twice, so all this sums up to –f(b)+2f(b) = +f(b), which is correct. Only the exit condition remains wrong.

Also n has to be even, not odd. The number of nodes between a and b is odd, since it equals n–1.

Code:

h = (b-a)/n
sum = f(a) + f(b)
k = 4

for i = 1 to n-1 do
  sum = sum + k * f(a + i*h)
  k = 6 - k
next

result = sum * h / 3

(12-12-2015 04:24 PM)Namir Wrote:  Many implementations use two summations--one for odd terms and one for even terms.

Right – this has a big advantage: If both sums are calculated separately, adding another iteration with 2·n is very easy since half of the function calls have already been evaluated. So you can calculate approximations for n=2, 4, 8, 16, ... until the results agree to the desired accuracy. This way you also get an estimate for the approximation error.

Dieter