inverse_of_bernoulli_numbers.sf 932 B

1234567891011121314151617181920212223242526272829303132333435
  1. #!/usr/bin/ruby
  2. # Daniel "Trizen" Șuteu
  3. # Date: 29 June 2017
  4. # Edit: 27 March 2023
  5. # https://github.com/trizen
  6. # Inverse of Bernoulli numbers, based on the inverse of the following asymptotic formula:
  7. # |Bn| ~ 2 / (2*pi)^n * n!
  8. # Using Stirling's approximation for n!, we have:
  9. # |Bn| ~ 2 / (2*pi)^n * sqrt(2*pi*n) * (n/e)^n
  10. # This gives us the following inverse formula:
  11. # n ~ lgrt((|Bn| / (4*pi))^(1/(2*pi*e))) * 2*pi*e - 1/2
  12. # Where `lgrt(n)` is defined as:
  13. # lgrt(x^x) = x
  14. define e = Num.e
  15. define τ = Num.tau
  16. func inverse_bernoulli_W (n) {
  17. (log(n/2) - log(τ)) / LambertW((log(n/2) - log(τ)) / (τ*e)) - 1/2
  18. }
  19. func inverse_bernoulli_lgrt (n) {
  20. lgrt((n/(2*τ))**(1/(e*τ))) * e * τ - 1/2
  21. }
  22. var x = abs(bernreal(1000000))
  23. say inverse_bernoulli_W(x) #=> 999999.999999996521295786570230337488233833193417
  24. say inverse_bernoulli_lgrt(x) #=> 999999.999999996521295786570230337488233833193417