float accuracy when calculating Fibonacci numbers

[Yuval Lifshitz]

Dear community,

I ran into the following issue when started playing with go. Had the following small program to calculate Fibonacci numbers using the closed form:

package "fib"

import "math"

var sqrt_5 = math.Sqrt(5.0)
var psi = -1.0/math.Phi

// had to add this to solve accuracy issue
func round(x float64) uint64 {
    return uint64(math.Floor(x + 0.5))
}

// find the nth fibonacci number based on Binet's formula
// see here: https://en.wikipedia.org/wiki/Fibonacci_number#Closed-form_expression
func fibonacci(n uint) uint64 {
    float_n := float64(n)
    return round((math.Pow(math.Phi, float_n) - math.Pow(psi, float_n))/sqrt_5)
}

When I first tried without doing the "rounding" I did not get whole numbers as the output (similar program in C++ gave whole numbers without rounding).

Following test is failing without the call to "round()":

package "fib"

import "testing"

var results = map[uint]uint64 {
    0: 0,
    2: 1,
    10: 55,
    14: 377,
    20: 6765,
    29: 514229,
    33: 3524578,
    59: 956722026041,
}

func Test_fibonacci(t *testing.T) {
    for arg, expected := range results {
        actual := fibonacci(arg)
        if actual != expected {
            t.Error("Expected", expected, "got", actual)
        }
    }
}

Are there known accuracy issues with floating points in go?

[Ian Lance Taylor]

No.

I suggest that you show us your C code. Thanks.

Ian

[Michael Jones]

Yuval,

There are fundamental issues here.

1. That equation (de Moivre, Binet) is from the algebra of ideal numbers. Numbers of infinite precision. Not the realm of computer arithmetic. It works fine with double precision (go: float64, c/c++: double) up to F(75) but must fail for F(76) due to the limited precision of 64-bit floating point and has nothing to do with language.

F(76) = 3416454622906707 but the best we can do in 64 bits is 3416454622906706 even with a Pow() function good to +/-1 least significant bit.

2. Another difference between algebra and computer arithmetic (well...actually about the floor function) is that one of your two power terms is not needed. Psi < 1 so Psi^N is pretty small, so small, that it never changes the value of the rounded result. So you can just evaluate:

return round(math.Pow(math.Phi, float_n) / sqrt_5)

3. Using a map in your test is not wrong, but it means that the tests will be executed in a random order, which makes the printed errors a little less clear.

celeste:fib mtj$ go test -v 

=== RUN   Test_fibonacci

--- FAIL: Test_fibonacci (0.00s)

fib_test.go:104: 82 : Expected 61305790721611591 got 61305790721611584

fib_test.go:104: 84 : Expected 160500643816367088 got 160500643816367040

fib_test.go:104: 81 : Expected 37889062373143906 got 37889062373143896

fib_test.go:104: 87 : Expected 679891637638612258 got 679891637638612096

fib_test.go:104: 88 : Expected 1100087778366101931 got 1100087778366101632

fib_test.go:104: 83 : Expected 99194853094755497 got 99194853094755488

fib_test.go:104: 77 : Expected 5527939700884757 got 5527939700884756

fib_test.go:104: 86 : Expected 420196140727489673 got 420196140727489600

fib_test.go:104: 90 : Expected 2880067194370816120 got 2880067194370815488

fib_test.go:104: 91 : Expected 4660046610375530309 got 4660046610375529472

fib_test.go:104: 92 : Expected 7540113804746346429 got 7540113804746344448

fib_test.go:104: 76 : Expected 3416454622906707 got 3416454622906706

fib_test.go:104: 85 : Expected 259695496911122585 got 259695496911122528

fib_test.go:104: 89 : Expected 1779979416004714189 got 1779979416004713728

fib_test.go:104: 79 : Expected 14472334024676221 got 14472334024676218

fib_test.go:104: 80 : Expected 23416728348467685 got 23416728348467676

FAIL

exit status 1

FAIL fib 0.003s

updated fib.go: https://play.golang.org/p/N-8lmjrMYAq

update fib_test.go: https://play.golang.org/p/FPuN58m1VVs

4. Plenty of ways to do this computation in Go using 64-bit integers, or big floats, or big ints.

5. Plenty of algorithms, some quite fascinating! https://ir.library.oregonstate.edu/downloads/t435gg51w

--
You received this message because you are subscribed to the Google Groups "golang-nuts" group.
To unsubscribe from this group and stop receiving emails from it, send an email to golang-nuts+unsubscribe@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.

--

Michael T. Jones
michae...@gmail.com

[Yuval Lifshitz]

(1) I understand the issue of limited precision, this is why I did not try anything above F(59) But my concern was not the difference between algebra and the go implementation it was the different results I got with the C/C++ implementation (gcc 7.3.1):

#include <math.h>

const double sqrt_5 = sqrt(5.0);
const double phi = (1.0 + sqrt_5)/2.0;
const double psi = -1.0/phi;

// find the nth fibonacci number based on Binet's formula
// see here: https://en.wikipedia.org/wiki/Fibonacci_number#Closed-form_expression

unsigned long long fibonacci(unsigned n)
{
    return (unsigned long long)((pow(phi, n) - pow(psi, n))/sqrt_5);
}

(2) Thanks for the fix, work well for the go implementation. However, if i try to apply it on the C++ one (where there is no rounding), it makes some of the tests fails. So, I guess that without rounding, the psi part is needed

(3) Is there a way to print the map ordered (didn't see that in your snippet)?

(4) Wanted to compare "double" in C++ to "float64" in go. As they should consume same amount of memory, have similar algorithms etc. but good to know about the "big" package in go

(5) Thanks for the reference. interesting read

To unsubscribe from this group and stop receiving emails from it, send an email to golang-nuts...@googlegroups.com.

For more options, visit https://groups.google.com/d/optout.

[Michael Jones]

yes. truncation and round-up are important differences between your c/c++ code and the go code.

[Yuval Lifshitz]

Thanks Michael. Is this "by design"? Does it worth an issue here?

[Louki Sumirniy]

I could be wrong but float64 is 'float' and 'double float' is float128. I dunno, I have not tinkered with these imprecise math types. That's just, essentially, what they were back in the olden days, like 20 years ago. 64 bit was basic float, 128 was double. I have been out of the loop for way too long on these things (and I am not a fan of floats anyway, give me some fixed precision pls).

[Yuval Lifshitz]

in C/C++ (at least on my system) "float" 4 has 4 bytes (32bits) and double has 8 bytes (64bits). according to this: https://golang.org/pkg/math/ its is the same thing with go.

[Michael Jones]

Not sure what you mean with "by design." The design here is the IEEE Standard for Floating-Point Arithmetic (IEEE 754). It is authoritative for Intel, AMD, IBM, Sun, etc.

Typical (universal?) C and C++ implementations map 'float' to 32-bit IEEE 754 binary floating point and 'double' to the related 64-bit version. Go uses this same hardware in the same (universal?) way, through the types "float32" and "float64."

Truncation of floating point types to integer, or the floor() function, discards any fractional part. Such that:

int(e) => 2

int(pi) => 3

Rounding up to the nearest integer with int(v+0.5) and floor(v+0.5) is a different operation:

round(e) => 3

round(pi) => 3

Your C program truncates and your Go program rounds. They are not the same program. They compute different things, which is why the results are different.

In a de Moivre coding of Fibonacci evaluation using rounding as above, Knuth showed that the pow((1-sqrt(5))/2,n) term--the second summand--may be dropped.

For 64-bit floating point the F(75) value is the most that you can generate this way.

Using other methods you can (and more quickly) generate values up to F(92), which is the greatest Fibonacci value (7540113804746346429) representable as a 64-bit unsigned integer.

For completeness, there is also the Ceiling function:

ceil(e) => 3

ceil(pi) => 4

hope this helps.

[Yuval Lifshitz]

It seems like go and C++ are doing something different regarding floating point arithmetic, but I cannot say one is better than the other. It is just the that C++ consistently overshoots (hence truncation work), and go consistently undershoots (this is why rounding is needed). For example, in C++:

F(59) = 956722026041.002

And in go, it is:

F(59) = 956722026040.999878

And for example:

[Marvin Renich]

* Yuval Lifshitz <yuv...@gmail.com> [180430 02:02]:

> It seems like go and C++ are doing something different regarding floating
> point arithmetic, but I cannot say one is better than the other. It is just
> the that C++ consistently overshoots (hence truncation work), and go
> consistently undershoots (this is why rounding is needed). For example, in
> C++:
>
> F(59) = 956722026041.002
>
> And in go, it is:
>
> F(59) = 956722026040.999878

/* aside - this doesn't seem relevant to this problem except perhaps
* the implementation of math.Pow, but...

The 8087 and successor FPUs (e.g. amd64 class processors) have an 80-bit
extended precision format which is used for all calculations. More than
a couple decades ago, when I was working on implementation of floating
point in an APL interpreter (written in C, C++, asm), the IEEE 754
standard did not specify how compilers should handle precision of
intermediate results, and the compiler we were using would only convert
from 80-bit to 64-bit when it needed to store a result back into memory.
This could easily give different results depending on how the floating
point expression was broken up. E.g. a = (b + c) - d could give a
different result than a = b + c; a = a - d.

The standard has, since then, attempted to address this, but I have not
studied this at all, so I don't know what is expected of C or Go.

*/

Back to your specific problem:

Go: psi = -0.61803398874989490
C++: psi = -0.61803398874989479

Go: math.Pow(math.Phi, float64(33)) = 7.88119600000012666e+06
C++: pow(phi, 33) = 7.88119600000013597e+06

And using hand-typed constants:

Go: math.Pow(1.61803398874989490, 33.0) = 7.88119600000012666e+06
C++: pow(1.61803398874989490, 33.0) = 7.88119600000013597e+06

So indeed there is a difference between math.Pow in Go and pow in C++.
I would say to file an issue. I have no idea which answer is closer to
correct, but using math/big and multiplying 1.61803398874989490 by
itself 33 times with a result precision of 80, I get
7.88119600000013562e+06.

Note: in Go, math.Phi is an untyped constant with precision greater
than float64. Adding var phi float64 = math.Phi and replacing
math.Phi everywhere else with phi did not change anything except the
value of psi (which then matched the C++ value). The C++ version of psi
is definitely farther from the "true" value, but as Michael said, the
error in pow(psi, n) is swamped by the value of pow(phi, n).

...Marvin

[andrey mirtchovski]

every time float accuracy is mentioned i think of this:

http://0.30000000000000004.com/

[Michael Jones]

Andrey, that's great!

On the Fibonacci series evaluation, let's make sure that we're all doing the same calculation. For completeness and safety, let's skip all library values and derived values. Here are more-than-sufficient versions of the three constants in Yuval's code:

// 40-decimal digit approximations

const sqrt5 float64 =  2.236067977499789696409173668731276235441  // math.Sqrt(5.0)

const phi   float64 =  1.618033988749894848204586834365638117720  // (1 + sqrt5) / 2

const psi   float64 = -0.6180339887498948482045868343656381177203 // (1 - sqrt5) / 2

Here is the evaluation: (https://play.golang.org/p/InwNEeyv4Bx)

func fibonacci(n uint) uint64 {

float_n := float64(n)

a := math.Pow(phi, float_n)

b := math.Pow(psi, float_n)

rounded := math.Floor(a/sqrt5 + 0.5)

both := (a - b) / sqrt5

fmt.Printf("%4d: rounded = %19.0f, both = %40.20f, delta = %+27.20e\n", n, rounded, both, both-rounded)

return uint64(rounded)

}

And the result:

   1: rounded =                   1, both =                   1.00000000000000000000, delta = +0.00000000000000000000e+00

   2: rounded =                   1, both =                   1.00000000000000000000, delta = +0.00000000000000000000e+00

   3: rounded =                   2, both =                   2.00000000000000000000, delta = +0.00000000000000000000e+00

   4: rounded =                   3, both =                   3.00000000000000000000, delta = +0.00000000000000000000e+00

   5: rounded =                   5, both =                   5.00000000000000000000, delta = +0.00000000000000000000e+00

   6: rounded =                   8, both =                   8.00000000000000000000, delta = +0.00000000000000000000e+00

   7: rounded =                  13, both =                  13.00000000000000000000, delta = +0.00000000000000000000e+00

   8: rounded =                  21, both =                  21.00000000000000000000, delta = +0.00000000000000000000e+00

   9: rounded =                  34, both =                  34.00000000000000000000, delta = +0.00000000000000000000e+00

  10: rounded =                  55, both =                  54.99999999999999289457, delta = -7.10542735760100185871e-15

  11: rounded =                  89, both =                  89.00000000000000000000, delta = +0.00000000000000000000e+00

  12: rounded =                 144, both =                 143.99999999999997157829, delta = -2.84217094304040074348e-14

  13: rounded =                 233, both =                 232.99999999999994315658, delta = -5.68434188608080148697e-14

  14: rounded =                 377, both =                 377.00000000000005684342, delta = +5.68434188608080148697e-14

  15: rounded =                 610, both =                 610.00000000000000000000, delta = +0.00000000000000000000e+00

  16: rounded =                 987, both =                 986.99999999999977262632, delta = -2.27373675443232059479e-13

  17: rounded =                1597, both =                1596.99999999999977262632, delta = -2.27373675443232059479e-13

  18: rounded =                2584, both =                2584.00000000000000000000, delta = +0.00000000000000000000e+00

  19: rounded =                4181, both =                4181.00000000000000000000, delta = +0.00000000000000000000e+00

  20: rounded =                6765, both =                6764.99999999999909050530, delta = -9.09494701772928237915e-13

  21: rounded =               10946, both =               10945.99999999999818101060, delta = -1.81898940354585647583e-12

  22: rounded =               17711, both =               17710.99999999999636202119, delta = -3.63797880709171295166e-12

  23: rounded =               28657, both =               28656.99999999999636202119, delta = -3.63797880709171295166e-12

  24: rounded =               46368, both =               46367.99999999999272404239, delta = -7.27595761418342590332e-12

  25: rounded =               75025, both =               75025.00000000000000000000, delta = +0.00000000000000000000e+00

  26: rounded =              121393, both =              121392.99999999998544808477, delta = -1.45519152283668518066e-11

  27: rounded =              196418, both =              196418.00000000000000000000, delta = +0.00000000000000000000e+00

  28: rounded =              317811, both =              317811.00000000000000000000, delta = +0.00000000000000000000e+00

  29: rounded =              514229, both =              514228.99999999994179233909, delta = -5.82076609134674072266e-11

  30: rounded =              832040, both =              832039.99999999988358467817, delta = -1.16415321826934814453e-10

  31: rounded =             1346269, both =             1346268.99999999976716935635, delta = -2.32830643653869628906e-10

  32: rounded =             2178309, both =             2178309.00000000000000000000, delta = +0.00000000000000000000e+00

  33: rounded =             3524578, both =             3524577.99999999953433871269, delta = -4.65661287307739257812e-10

  34: rounded =             5702887, both =             5702886.99999999906867742538, delta = -9.31322574615478515625e-10

  35: rounded =             9227465, both =             9227465.00000000000000000000, delta = +0.00000000000000000000e+00

  36: rounded =            14930352, both =            14930351.99999999813735485077, delta = -1.86264514923095703125e-09

  37: rounded =            24157817, both =            24157816.99999999627470970154, delta = -3.72529029846191406250e-09

  38: rounded =            39088169, both =            39088168.99999999254941940308, delta = -7.45058059692382812500e-09

  39: rounded =            63245986, both =            63245985.99999999254941940308, delta = -7.45058059692382812500e-09

  40: rounded =           102334155, both =           102334154.99999998509883880615, delta = -1.49011611938476562500e-08

  41: rounded =           165580141, both =           165580140.99999997019767761230, delta = -2.98023223876953125000e-08

  42: rounded =           267914296, both =           267914295.99999994039535522461, delta = -5.96046447753906250000e-08

  43: rounded =           433494437, both =           433494437.00000000000000000000, delta = +0.00000000000000000000e+00

  44: rounded =           701408733, both =           701408732.99999988079071044922, delta = -1.19209289550781250000e-07

  45: rounded =          1134903170, both =          1134903169.99999976158142089844, delta = -2.38418579101562500000e-07

  46: rounded =          1836311903, both =          1836311903.00000000000000000000, delta = +0.00000000000000000000e+00

  47: rounded =          2971215073, both =          2971215072.99999952316284179688, delta = -4.76837158203125000000e-07

  48: rounded =          4807526976, both =          4807526975.99999904632568359375, delta = -9.53674316406250000000e-07

  49: rounded =          7778742049, both =          7778742048.99999809265136718750, delta = -1.90734863281250000000e-06

  50: rounded =         12586269025, both =         12586269024.99999809265136718750, delta = -1.90734863281250000000e-06

  51: rounded =         20365011074, both =         20365011074.00000000000000000000, delta = +0.00000000000000000000e+00

  52: rounded =         32951280099, both =         32951280098.99999237060546875000, delta = -7.62939453125000000000e-06

  53: rounded =         53316291173, both =         53316291172.99998474121093750000, delta = -1.52587890625000000000e-05

  54: rounded =         86267571272, both =         86267571271.99998474121093750000, delta = -1.52587890625000000000e-05

  55: rounded =        139583862445, both =        139583862444.99996948242187500000, delta = -3.05175781250000000000e-05

  56: rounded =        225851433717, both =        225851433716.99996948242187500000, delta = -3.05175781250000000000e-05

  57: rounded =        365435296162, both =        365435296161.99993896484375000000, delta = -6.10351562500000000000e-05

  58: rounded =        591286729879, both =        591286729878.99987792968750000000, delta = -1.22070312500000000000e-04

  59: rounded =        956722026041, both =        956722026040.99987792968750000000, delta = -1.22070312500000000000e-04

  60: rounded =       1548008755920, both =       1548008755919.99975585937500000000, delta = -2.44140625000000000000e-04

  61: rounded =       2504730781961, both =       2504730781960.99951171875000000000, delta = -4.88281250000000000000e-04

  62: rounded =       4052739537881, both =       4052739537880.99902343750000000000, delta = -9.76562500000000000000e-04

  63: rounded =       6557470319842, both =       6557470319841.99902343750000000000, delta = -9.76562500000000000000e-04

  64: rounded =      10610209857723, both =      10610209857722.99804687500000000000, delta = -1.95312500000000000000e-03

  65: rounded =      17167680177565, both =      17167680177564.99609375000000000000, delta = -3.90625000000000000000e-03

  66: rounded =      27777890035288, both =      27777890035287.99609375000000000000, delta = -3.90625000000000000000e-03

  67: rounded =      44945570212853, both =      44945570212852.99218750000000000000, delta = -7.81250000000000000000e-03

  68: rounded =      72723460248141, both =      72723460248140.98437500000000000000, delta = -1.56250000000000000000e-02

  69: rounded =     117669030460994, both =     117669030460993.98437500000000000000, delta = -1.56250000000000000000e-02

  70: rounded =     190392490709135, both =     190392490709134.96875000000000000000, delta = -3.12500000000000000000e-02

  71: rounded =     308061521170129, both =     308061521170129.00000000000000000000, delta = +0.00000000000000000000e+00

  72: rounded =     498454011879264, both =     498454011879263.87500000000000000000, delta = -1.25000000000000000000e-01

  73: rounded =     806515533049393, both =     806515533049392.87500000000000000000, delta = -1.25000000000000000000e-01

  74: rounded =    1304969544928657, both =    1304969544928656.75000000000000000000, delta = -2.50000000000000000000e-01

  75: rounded =    2111485077978050, both =    2111485077978049.50000000000000000000, delta = -5.00000000000000000000e-01

  76: rounded =    3416454622906706, both =    3416454622906706.00000000000000000000, delta = +0.00000000000000000000e+00

  77: rounded =    5527939700884756, both =    5527939700884755.00000000000000000000, delta = -1.00000000000000000000e+00

  78: rounded =    8944394323791464, both =    8944394323791463.00000000000000000000, delta = -1.00000000000000000000e+00

  79: rounded =   14472334024676218, both =   14472334024676218.00000000000000000000, delta = +0.00000000000000000000e+00

  80: rounded =   23416728348467676, both =   23416728348467676.00000000000000000000, delta = +0.00000000000000000000e+00

  81: rounded =   37889062373143896, both =   37889062373143896.00000000000000000000, delta = +0.00000000000000000000e+00

  82: rounded =   61305790721611584, both =   61305790721611584.00000000000000000000, delta = +0.00000000000000000000e+00

  83: rounded =   99194853094755488, both =   99194853094755488.00000000000000000000, delta = +0.00000000000000000000e+00

  84: rounded =  160500643816367040, both =  160500643816367040.00000000000000000000, delta = +0.00000000000000000000e+00

  85: rounded =  259695496911122528, both =  259695496911122528.00000000000000000000, delta = +0.00000000000000000000e+00

  86: rounded =  420196140727489600, both =  420196140727489600.00000000000000000000, delta = +0.00000000000000000000e+00

  87: rounded =  679891637638612096, both =  679891637638612096.00000000000000000000, delta = +0.00000000000000000000e+00

  88: rounded = 1100087778366101632, both = 1100087778366101632.00000000000000000000, delta = +0.00000000000000000000e+00

  89: rounded = 1779979416004713728, both = 1779979416004713728.00000000000000000000, delta = +0.00000000000000000000e+00

  90: rounded = 2880067194370815488, both = 2880067194370815488.00000000000000000000, delta = +0.00000000000000000000e+00

  91: rounded = 4660046610375529472, both = 4660046610375529472.00000000000000000000, delta = +0.00000000000000000000e+00

  92: rounded = 7540113804746344448, both = 7540113804746344448.00000000000000000000, delta = +0.00000000000000000000e+00

Yuval, please do the same test on your computer. Please evaluate the two powers as separate expressions shown here ("a" and "b"). Let's see if things agree.

On Mon, Apr 30, 2018 at 10:18 AM andrey mirtchovski <mirtc...@gmail.com> wrote:

every time float accuracy is mentioned i think of this:

http://0.30000000000000004.com/

--
You received this message because you are subscribed to the Google Groups "golang-nuts" group.
To unsubscribe from this group and stop receiving emails from it, send an email to golang-nuts...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.

[Marvin Renich]

* Michael Jones <michae...@gmail.com> [180430 13:54]:

> Andrey, that's great!
>
> On the Fibonacci series evaluation, let's make sure that we're all doing
> the same calculation. For completeness and safety, let's skip all library
> values and derived values. Here are more-than-sufficient versions of the
> three constants in Yuval's code:

I think both of you missed the point of my message. This has nothing to
do specifically with Fibonacci numbers. The problem is simply that
math.Pow does not give as precise answers as the C++ std lib. Look at
https://play.golang.org/p/_gVbAWjeoyW and notice that the result from
math.Pow is off in the 15th decimal digit. Both the bigPow and the C++
pow library function are good until the 17th digit, which is what is
expected. (Actually, I thought bigPow might do better, with the
specified precision of 300.) I trust the Mathematica answer to be
correct to the number of digits printed (last digit rounded).

There is no math.Phi involved here, except that the hand-typed constant
happens to be as close as you can get to Phi in a float64 (but this is
not relevant to the test).

...Marvin

[andrey mirtchovski]

> math.Pow does not give as precise answers as the C++ std lib.

math pow doesn't just multiply X by itself Y times, it uses a
different algorithm:

https://golang.org/src/math/pow.go#L40

it would perhaps make sense to quantify the error margins that this
algorithm introduces.

not sure what c++ stdlib does.

[Michael Jones]

using an optimally-precise Pow:

   1: rounded =                   1, both =                   1.00000000000000000000, delta = +0.00000000000000000000e+00

   2: rounded =                   1, both =                   1.00000000000000000000, delta = +0.00000000000000000000e+00

   3: rounded =                   2, both =                   2.00000000000000000000, delta = +0.00000000000000000000e+00

   4: rounded =                   3, both =                   3.00000000000000000000, delta = +0.00000000000000000000e+00

   5: rounded =                   5, both =                   5.00000000000000000000, delta = +0.00000000000000000000e+00

   6: rounded =                   8, both =                   8.00000000000000000000, delta = +0.00000000000000000000e+00

   7: rounded =                  13, both =                  13.00000000000000000000, delta = +0.00000000000000000000e+00

   8: rounded =                  21, both =                  21.00000000000000000000, delta = +0.00000000000000000000e+00

   9: rounded =                  34, both =                  34.00000000000000000000, delta = +0.00000000000000000000e+00

  10: rounded =                  55, both =                  55.00000000000000000000, delta = +0.00000000000000000000e+00

  11: rounded =                  89, both =                  89.00000000000000000000, delta = +0.00000000000000000000e+00

  12: rounded =                 144, both =                 144.00000000000000000000, delta = +0.00000000000000000000e+00

  13: rounded =                 233, both =                 233.00000000000000000000, delta = +0.00000000000000000000e+00

  14: rounded =                 377, both =                 377.00000000000005684342, delta = +5.68434188608080148697e-14

  15: rounded =                 610, both =                 610.00000000000000000000, delta = +0.00000000000000000000e+00

  16: rounded =                 987, both =                 986.99999999999977262632, delta = -2.27373675443232059479e-13

  17: rounded =                1597, both =                1596.99999999999977262632, delta = -2.27373675443232059479e-13

  18: rounded =                2584, both =                2584.00000000000000000000, delta = +0.00000000000000000000e+00

  19: rounded =                4181, both =                4181.00000000000000000000, delta = +0.00000000000000000000e+00

  20: rounded =                6765, both =                6765.00000000000090949470, delta = +9.09494701772928237915e-13

  21: rounded =               10946, both =               10946.00000000000363797881, delta = +3.63797880709171295166e-12

  22: rounded =               17711, both =               17711.00000000000000000000, delta = +0.00000000000000000000e+00

  23: rounded =               28657, both =               28657.00000000000363797881, delta = +3.63797880709171295166e-12

  24: rounded =               46368, both =               46368.00000000000727595761, delta = +7.27595761418342590332e-12

  25: rounded =               75025, both =               75025.00000000001455191523, delta = +1.45519152283668518066e-11

  26: rounded =              121393, both =              121393.00000000001455191523, delta = +1.45519152283668518066e-11

  27: rounded =              196418, both =              196418.00000000005820766091, delta = +5.82076609134674072266e-11

  28: rounded =              317811, both =              317811.00000000011641532183, delta = +1.16415321826934814453e-10

  29: rounded =              514229, both =              514228.99999999994179233909, delta = -5.82076609134674072266e-11

  30: rounded =              832040, both =              832040.00000000000000000000, delta = +0.00000000000000000000e+00

  31: rounded =             1346269, both =             1346269.00000000046566128731, delta = +4.65661287307739257812e-10

  32: rounded =             2178309, both =             2178309.00000000000000000000, delta = +0.00000000000000000000e+00

  33: rounded =             3524578, both =             3524577.99999999906867742538, delta = -9.31322574615478515625e-10

  34: rounded =             5702887, both =             5702886.99999999906867742538, delta = -9.31322574615478515625e-10

  35: rounded =             9227465, both =             9227465.00000000000000000000, delta = +0.00000000000000000000e+00

  36: rounded =            14930352, both =            14930352.00000000000000000000, delta = +0.00000000000000000000e+00

  37: rounded =            24157817, both =            24157816.99999999254941940308, delta = -7.45058059692382812500e-09

  38: rounded =            39088169, both =            39088169.00000000000000000000, delta = +0.00000000000000000000e+00

  39: rounded =            63245986, both =            63245986.00000000745058059692, delta = +7.45058059692382812500e-09

  40: rounded =           102334155, both =           102334155.00000002980232238770, delta = +2.98023223876953125000e-08

  41: rounded =           165580141, both =           165580141.00000002980232238770, delta = +2.98023223876953125000e-08

  42: rounded =           267914296, both =           267914296.00000014901161193848, delta = +1.49011611938476562500e-07

  43: rounded =           433494437, both =           433494437.00000011920928955078, delta = +1.19209289550781250000e-07

  44: rounded =           701408733, both =           701408733.00000023841857910156, delta = +2.38418579101562500000e-07

  45: rounded =          1134903170, both =          1134903170.00000047683715820312, delta = +4.76837158203125000000e-07

  46: rounded =          1836311903, both =          1836311903.00000071525573730469, delta = +7.15255737304687500000e-07

  47: rounded =          2971215073, both =          2971215072.99999904632568359375, delta = -9.53674316406250000000e-07

  48: rounded =          4807526976, both =          4807526976.00000190734863281250, delta = +1.90734863281250000000e-06

  49: rounded =          7778742049, both =          7778742049.00000286102294921875, delta = +2.86102294921875000000e-06

  50: rounded =         12586269025, both =         12586269025.00000572204589843750, delta = +5.72204589843750000000e-06

  51: rounded =         20365011074, both =         20365011074.00000762939453125000, delta = +7.62939453125000000000e-06

  52: rounded =         32951280099, both =         32951280099.00001144409179687500, delta = +1.14440917968750000000e-05

  53: rounded =         53316291173, both =         53316291173.00003051757812500000, delta = +3.05175781250000000000e-05

  54: rounded =         86267571272, both =         86267571272.00003051757812500000, delta = +3.05175781250000000000e-05

  55: rounded =        139583862445, both =        139583862445.00006103515625000000, delta = +6.10351562500000000000e-05

  56: rounded =        225851433717, both =        225851433717.00018310546875000000, delta = +1.83105468750000000000e-04

  57: rounded =        365435296162, both =        365435296162.00024414062500000000, delta = +2.44140625000000000000e-04

  58: rounded =        591286729879, both =        591286729878.99987792968750000000, delta = -1.22070312500000000000e-04

  59: rounded =        956722026041, both =        956722026041.00073242187500000000, delta = +7.32421875000000000000e-04

  60: rounded =       1548008755920, both =       1548008755920.00024414062500000000, delta = +2.44140625000000000000e-04

  61: rounded =       2504730781961, both =       2504730781960.99902343750000000000, delta = -9.76562500000000000000e-04

  62: rounded =       4052739537881, both =       4052739537881.00341796875000000000, delta = +3.41796875000000000000e-03

  63: rounded =       6557470319842, both =       6557470319842.00585937500000000000, delta = +5.85937500000000000000e-03

  64: rounded =      10610209857723, both =      10610209857722.99804687500000000000, delta = -1.95312500000000000000e-03

  65: rounded =      17167680177565, both =      17167680177564.99609375000000000000, delta = -3.90625000000000000000e-03

  66: rounded =      27777890035288, both =      27777890035287.99218750000000000000, delta = -7.81250000000000000000e-03

  67: rounded =      44945570212853, both =      44945570212852.99218750000000000000, delta = -7.81250000000000000000e-03

  68: rounded =      72723460248141, both =      72723460248140.98437500000000000000, delta = -1.56250000000000000000e-02

  69: rounded =     117669030460994, both =     117669030460993.98437500000000000000, delta = -1.56250000000000000000e-02

  70: rounded =     190392490709135, both =     190392490709135.00000000000000000000, delta = +0.00000000000000000000e+00

  71: rounded =     308061521170129, both =     308061521170129.18750000000000000000, delta = +1.87500000000000000000e-01

  72: rounded =     498454011879264, both =     498454011879264.00000000000000000000, delta = +0.00000000000000000000e+00

  73: rounded =     806515533049393, both =     806515533049393.12500000000000000000, delta = +1.25000000000000000000e-01

  74: rounded =    1304969544928657, both =    1304969544928656.50000000000000000000, delta = -5.00000000000000000000e-01

  75: rounded =    2111485077978050, both =    2111485077978050.00000000000000000000, delta = +0.00000000000000000000e+00

  76: rounded =    3416454622906708, both =    3416454622906707.50000000000000000000, delta = -5.00000000000000000000e-01

  77: rounded =    5527939700884756, both =    5527939700884756.00000000000000000000, delta = +0.00000000000000000000e+00

  78: rounded =    8944394323791466, both =    8944394323791466.00000000000000000000, delta = +0.00000000000000000000e+00

  79: rounded =   14472334024676214, both =   14472334024676214.00000000000000000000, delta = +0.00000000000000000000e+00

  80: rounded =   23416728348467700, both =   23416728348467700.00000000000000000000, delta = +0.00000000000000000000e+00

  81: rounded =   37889062373143904, both =   37889062373143904.00000000000000000000, delta = +0.00000000000000000000e+00

  82: rounded =   61305790721611608, both =   61305790721611608.00000000000000000000, delta = +0.00000000000000000000e+00

  83: rounded =   99194853094755536, both =   99194853094755536.00000000000000000000, delta = +0.00000000000000000000e+00

  84: rounded =  160500643816367264, both =  160500643816367264.00000000000000000000, delta = +0.00000000000000000000e+00

  85: rounded =  259695496911122752, both =  259695496911122752.00000000000000000000, delta = +0.00000000000000000000e+00

  86: rounded =  420196140727489984, both =  420196140727489984.00000000000000000000, delta = +0.00000000000000000000e+00

  87: rounded =  679891637638612608, both =  679891637638612608.00000000000000000000, delta = +0.00000000000000000000e+00

  88: rounded = 1100087778366102528, both = 1100087778366102528.00000000000000000000, delta = +0.00000000000000000000e+00

  89: rounded = 1779979416004716032, both = 1779979416004716032.00000000000000000000, delta = +0.00000000000000000000e+00

  90: rounded = 2880067194370818048, both = 2880067194370818048.00000000000000000000, delta = +0.00000000000000000000e+00

  91: rounded = 4660046610375532544, both = 4660046610375532544.00000000000000000000, delta = +0.00000000000000000000e+00

  92: rounded = 7540113804746352640, both = 7540113804746352640.00000000000000000000, delta = +0.00000000000000000000e+00

--
You received this message because you are subscribed to the Google Groups "golang-nuts" group.
To unsubscribe from this group and stop receiving emails from it, send an email to golang-nuts...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.

[Marvin Renich]

* andrey mirtchovski <mirtc...@gmail.com> [180430 17:41]:

> > math.Pow does not give as precise answers as the C++ std lib.
>
> math pow doesn't just multiply X by itself Y times, it uses a
> different algorithm:
>
> https://golang.org/src/math/pow.go#L40

Sure, I understand that. However, the current algorithm used by
math.Pow produces results that are far enough from the true value that
it is worth opening an issue.

> it would perhaps make sense to quantify the error margins that this
> algorithm introduces.

As far as the OP's question about why he needs to round in Go but
truncate in C++, I am not convinced that the error in math.Pow is
significant.

Binet's formula, given exact (non-computer limited) values and
calculations, will give exact integral values for the Fibonacci numbers.
If you start with the assumption that using Binet's formula with IEEE
floating point arithmetic will give results that are close enough to the
correct value that you will have no trouble determining which integer is
the right one, then the assumption that the error is between -0.5 and
+0.5 is much more reasonable than the assumption that the error is in
the range [0, 1).

It is certainly plausible that one could prove that the error will
always be positive, given known characteristics of the computer
program's implementation of pow, among other things, but without such
analysis, assuming that truncating the result will always give the
correct answer is naive.

I posit that the OP's C++ implementation only works because of specific
characteristics of C++ library pow function, the last-mantissa-bit
rounding of the value used for phi, and other aspects of the
implementation that have not be adequately analyzed to be certain that
the error is always positive.

In short, the OP's assumption that the algorithm in his C++ program is
correct is wrong; it just happens to work for the values tested. Both
the C++ and the Go program should round the floating point result to get
the correct integer.

...Marvin

[Paul Hankin]

On Friday, 27 April 2018 23:57:42 UTC+2, Michael Jones wrote:

Yuval,

There are fundamental issues here.

1. That equation (de Moivre, Binet) is from the algebra of ideal numbers. Numbers of infinite precision. Not the realm of computer arithmetic. It works fine with double precision (go: float64, c/c++: double) up to F(75) but must fail for F(76) due to the limited precision of 64-bit floating point and has nothing to do with language.

F(76) = 3416454622906707 but the best we can do in 64 bits is 3416454622906706 even with a Pow() function good to +/-1 least significant bit.

2. Another difference between algebra and computer arithmetic (well...actually about the floor function) is that one of your two power terms is not needed. Psi < 1 so Psi^N is pretty small, so small, that it never changes the value of the rounded result. So you can just evaluate:

return round(math.Pow(math.Phi, float_n) / sqrt_5)

A neat not-very-well-known trick is to use that phi^n = Fib(n)phi + Fib(n-1). Numbers of the form a*phi+b where a and b are integers are closed under multiplication, so you can compute phi^n exactly using integer arithmetic -- in log(n) arithmetic operations.

https://play.golang.org/p/j4mZ93c820R

-- 

Paul

[Yuval Lifshitz]

seems like math accuracy issues are known, or even, by design. see:

https://github.com/golang/go/issues/9546

https://github.com/golang/go/issues/9545