Post Reply 
(49g 50g) Some Relatives of Bernoulli Numbers
05-09-2022, 03:25 PM (This post was last modified: 05-09-2022 03:28 PM by John Keith.)
Post: #1
(49g 50g) Some Relatives of Bernoulli Numbers
This is a collection of three programs which compute Poly-Bernoulli numbers, Genocchi numbers and Gregory coefficients, numbers related to Bernoulli numbers.

Given an integer n on level 2 and a number k on level 1, the first program PBERN returns the Poly-Bernoulli number Bn(k). If k is a negative integer, the result will be an integer. Otherwise it will be a rational number (a fraction). The level 1 argument may be symbolic, but the resulting expression may be very large.

Given an integer n, the second program GENOC will return n terms of the unsigned Genocchi numbers to level 2 and the Genocchi medians A005439 to level 1. The reason for the program returning both sequences is that the Seidel algorithm used by the program generates both at once.

The last program Gn returns a list of the Gregory coefficients from 0 to n. Also known as the Bernoulli numbers of the second kind, they are rational numbers as well.

All three programs require the ListExt library. The programs are posted as a directory for convenience but each program is independent and any unused programs may be deleted or moved to other directories.

Code:

DIR
  PBERN
  \<< OVER I\->R \-> n k r
    \<< r 1.
      IF >
      THEN { 1 } 2. r
        START 0 SWAP + 2.
          \<< NSUB R\->I * +
          \>> DOSUBS 1 +          @ Stirling S2 numbers
        NEXT n LSEQ :: * LSCAN
        -1 -1 n LMSEQ * *         @ Multiply by signed factorials
        2 n 1 + LSEQR k ^ 2.      @ {2..n}^k
        \<< / EVAL
        \>> DOLIST \GSLIST EVAL   @ Divide, sum and simplify
        r 2. MOD :: NEG IFT EVAL  @ Negate if n is odd
      ELSE r 1. SAME
        { 2 k ^ INV EVAL }        @ Return 1/2^k if n = 1
        { n 1 + } IFTE            @ Return 1 if n = 0, preserving type
      END
    \>>
  \>>
  GENOC
  \<< I\->R \-> n
    \<< 1 1 { 1 } 3. n 2. *
      FOR k
        IF k 2. MOD
        THEN 0 + :: + LSCAN DUP LPOPR NIP
        ELSE REV :: + LSCAN REV DUP HEAD      @ Seidel algorithm
        END NEWOB SWAP
      NEXT DROP n 2. * \->LIST -2 LDIST EVAL  @ Break into 2 lists
    \>>
  \>>
  Gn
  \<< 1 - 2 OVER 1 + LSEQR :: * LSCAN \-> n f  @ Factorials 2..n
    \<< 1 2 INV { 1 } 1 n
      FOR k 0 SWAP + 2.
        \<< k * -
        \>> DOSUBS 1 + DUP 2 k 2 + LSEQR /     @ Stirling S1 numbers
          \GSLIST f k GET / EVAL SWAP
      NEXT DROP n 2 + \->LIST
    \>>
  \>>
END
Find all posts by this user
Quote this message in a reply
Post Reply 




User(s) browsing this thread: 1 Guest(s)