HP 50g Romberg Integration

[Gerson W. Barbosa]

(04-02-2014 10:02 PM)C.Ret Wrote:  I just recover from an old notepad the ROMBERG INTEGRATION programm I wrote for my HP-28S a few ago ...

It's regular RPL code, so, it will also run on any "new" RPL systems such as the HP-50g !

As you said, your code runs with no modification on the HP 50g:

Code:


%%HP: T(3)A(D)F(.);
\<< \-> Fct a b m
  \<< { m m } 0. CON 1. m
    FOR i 1. i
      FOR j
        IF j 1. ==
        THEN
          IF i 1. ==
          THEN b a - 2. / a Fct EVAL b Fct EVAL + *
          ELSE DUP i 1. - { 1. } + GET 2. / b a - 2. i 1. - ^ / 0. 1. 2. i 1. - ^
            FOR k OVER k * a + Fct EVAL + 2.
            STEP * +
          END
        ELSE DUP i 1. - j 1. - 2. \->LIST GET OVER i j 1. - 2. \->LIST GET 4. j 1. - ^ SWAP OVER * ROT - SWAP 1. - /
        END { i j } SWAP PUT
      NEXT
    NEXT { m m } GET
  \>>
\>>

A test:

Code:


« LN NEG LN NEG »
0.000000000001
0.999999999999
20
« ROMBERG » TEVAL

0.577222197843
s:290.1873               (4 m 50.1 s on the emulator --> about three hours on the HP 50g)

The buit-in integrator, OTOH:

Code:


« '∫(0.000000000001,0.999999999999,-LN(-LN(X)),X)' EVAL » TEVAL
0.577215664728
s: 52.8994               (52.9 s on the emulator)


« '∫(0.,1.,-LN(-LN(X)),X)' EVAL » TEVAL
0.577215664748
s:1934.7719              (32 m 15.5 s on the HP 50g)

The actual result is 0.577215664902 ( « 1. Psi NEG» ). I understand your program and peacecalc's are not meant to replace the built-in integrator, but I've read the legacy algorithm is not the best one around (it's only there because it's been tried and tested since the HP-34C, the Saturn code since the HP-71B, I think). Does anyone know of alternative programs that have better performance on this example, that is, faster and more accurate? Thanks!

Gerson.