inverse_of_bernoulli_numbers.sf

#!/usr/bin/ruby

# Daniel "Trizen" Șuteu
# Date: 29 June 2017
# Edit: 27 March 2023
# https://github.com/trizen

# Inverse of Bernoulli numbers, based on the inverse of the following asymptotic formula:
#   |Bn| ~ 2 / (2*pi)^n * n!

# Using Stirling's approximation for n!, we have:
#   |Bn| ~ 2 / (2*pi)^n * sqrt(2*pi*n) * (n/e)^n

# This gives us the following inverse formula:
#   n ~ lgrt((|Bn| / (4*pi))^(1/(2*pi*e))) * 2*pi*e - 1/2

# Where `lgrt(n)` is defined as:
#   lgrt(x^x) = x

define e = Num.e
define τ = Num.tau

func inverse_bernoulli_W (n) {
    (log(n/2) - log(τ)) / LambertW((log(n/2) - log(τ)) / (τ*e)) - 1/2
}

func inverse_bernoulli_lgrt (n) {
    lgrt((n/(2*τ))**(1/(e*τ))) * e * τ - 1/2
}

var x = abs(bernreal(1000000))

say inverse_bernoulli_W(x)         #=> 999999.999999996521295786570230337488233833193417
say inverse_bernoulli_lgrt(x)      #=> 999999.999999996521295786570230337488233833193417