Bernoulli function in Haskell

Question

I want to write a bernoulli function bernoulli:: Integer -> Rational in haskell, using the following algorithm for calculating the bernoulli number for a given integer.

the functions "frac" and "binom" are used to calculate the binomial in the definition. this is what i have so far:

fact :: Integer -> Integer
fact i = foldr (*) 1 [1..i]

binom :: Integer -> Integer -> Integer
binom n k = (fact n) `div` (fact k* fact (n-k))

bernoulli :: Integer -> Rational
bernoulli 0 = 1
bernoulli i = ((binom i j) * (bernoulli j)) / (j - i - 1) where j = i-1

I've tried it a few different times now, but either the recursion doesn't work or the resulting Rational is wrong.

Can you post the wire associates with the code? Right now it reads as "something's wrong" which is vague. — Thomas M. DuBuisson, Dec 02 '21 at 01:44
Note that this definition of `binom` is horribly inefficient. Most of the terms in the fraction cancel. It's almost never a good idea to actually compute factorials, instead take the product over only the terms that don't cancel. — leftaroundabout, Dec 02 '21 at 08:10

David Lukas · Accepted Answer · 2021-12-02T11:30:23.293

I found three problems in your code:

Parentheses in the binom
Mixing Rational and Integer
Your function bernoulli is not a sum, but only first member

In my code, you can see how I coped with these problems.

fact :: Integer -> Integer
fact i = foldr (*) 1 [1..i]

binom :: Integer -> Integer -> Integer
binom n k = (fact n) `div` ((fact k) * fact (n-k))

bernoulli :: Integer -> Rational
bernoulli 0 = 1
bernoulli n = sum [
    toRational(binom n k) * (bernoulli k) / toRational(k - n - 1)
    | k <- [0..(n-1)]
  ]

Test:

map bernoulli [0..10]

Output:

[1 % 1,(-1) % 2,1 % 6,0 % 1,(-1) % 30,0 % 1,1 % 42,0 % 1,(-1) % 30,0 % 1,5 % 66]

Small addition:
If we don't follow the rule of leveraging existing libraries, the solution could also look like this:

binom :: Rational -> Rational -> Rational
binom n k = product [ ( n + 1 - i ) / i | i <- [ 1 .. k ] ]

bernoulli :: Rational -> Rational
bernoulli 0 = 1
bernoulli n = sum [
    binom n k * bernoulli k / (k - n - 1)
    | k <- [0..(n-1)]
  ]

Because:

Note the similarity of the program with mathematical notation.

p.s. bernoulli with foldr:

bernoulli :: Integer -> Rational
bernoulli n = foldr (summand n) (0::Rational) [0 .. (n-1)]
  where
    summand n k s1 = s1 + toRational(binom n k) * bernoulli k / toRational(k - n - 1)

The parentheses are mostly superfluous. – leftaroundabout Dec 02 '21 at 08:12 — leftaroundabout, Dec 02 '21 at 08:12

Bernoulli function in Haskell

1 Answers1