Simpson's 1/3 rule on the 50g

[EngineerX]

Hello,

I would like to implement the following on the 50g in RPL but don't know how. Seems simple but I am new to the language. I have a feeling the Prime would be a better choice for me who has some experience with Matlab, but the thing is that I prefer the 50g for nostalgic reasons, although I never really learned more than the trivial.

I would like to be able to calculate the area under a curve using Simpson's 1/3 rule for a given set of data.

For exemple:

For x = {x0 x1 x2 x3 ... xN} where N is an even number.

Compute
A = h/3*(y0 + 4*y1 + 2*y2 + 4*y3 + 2*y4 + 4*y5 +...+ yN)

Very much appreciated.

[Wes Loewer]

I dug through my old files an found a program that I wrote some 15+ years ago. Looking at it now, I might have done this differently, but here it is, warts and all.

To understand the code fully, look up the syntax of the Σ and the | (where) commands.

Code:


@ SIMPSONS(xleft,xright,fnc,xvar,steps)
@ calculates area via Simpson's Rule (using parabolas)
@ Wes Loewer
@
@ parameters:
@ xleft   - left bound
@ xright  - right bound
@ fnc     - function
@ xvar   - independent variable
@ steps   - number of steps, must be even (ie, mult. of 2)

« → xleft xright fnc xvar steps
  «
    PUSH                      @ save flags
    -105. SF                  @ approx mode ON
    xleft →NUM 'xleft' STO    @ convert to decimals
    xright →NUM 'xright' STO
                              @ make sure steps is an even integer
    steps →NUM IP 2. / CEIL 2. * 'steps' STO

    xright xleft - steps / 
    → deltax                  @ step width
    «
      @ Simpson's Rule is
      @ (d/3)*( f(x0) + 4(f(x1) + 2f(x2) + 4f(x3) + 2f(x4) + ...) + f(xn) )

      @ middle odd terms * 4
      'n' 1. steps 2. /       @ 1st 3 arg of Σ_{n=1}^{steps/2}
                              @ f(xleft + (2n-1)*deltax)
      fnc xvar xleft '2.*n-1.' deltax * + 2. →LIST |
      Σ 4. *                  @ sum up middle odd terms *4

      @ middle even terms * 2
      'n' 1. steps 2. / 1. -  @ 1st 3 arg of Σ_{n=1}^{steps/2-1}
                              @ f(xleft + (2n)*deltax)
      fnc xvar xleft '2.*n' deltax * + 2. →LIST |
      @ Σ correctly handles no terms, Σ_{n=1}^{0}
      Σ 2. *                  @ sum up middle even terms *2
      +

      @ end terms
      fnc xvar xleft  2. →LIST | →NUM @ f(xleft)
      +
      fnc xvar xright 2. →LIST | →NUM @ f(xright)
      +

      deltax * 3. /           @ * deltax / 3 to get area
    »
    POP                       @ restore flags
  »
»

Feel free to ask if there is anything that is not clear.

[Juan14]

This program uses GETI to get trough each element of the list, flag -64 is checked in the WHILE loop, it gets set when the index reach the end of the list. In the loop, a 4 or 2 is generated with the help of MOD. The program requires h in level 2: and the list of values in level 1:

«
1 GETI UNROT
WHILE GETI -64 FC?
REPEAT
OVER 2 MOD 2 * 2 + * 4 ROLL + UNROT
END NIP NIP + * 3 /
»

[EngineerX]

Thank you for the inputs. Simple and elegant. I would have used many more commands in attempt to emulate what Matlab does.