[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
[Octave-bug-tracker] [bug #60682] betainc is inaccurate
From: |
Nicholas Jankowski |
Subject: |
[Octave-bug-tracker] [bug #60682] betainc is inaccurate |
Date: |
Wed, 11 Aug 2021 19:47:07 -0400 (EDT) |
User-agent: |
Mozilla/5.0 (Windows NT 10.0; Win64; x64) AppleWebKit/537.36 (KHTML, like Gecko) Chrome/92.0.4515.131 Safari/537.36 |
Follow-up Comment #4, bug #60682 (project octave):
curious how matlab compares, here are the outputs of the two commands below in
Matlab 2021a:
>> x=2.^(-(0:20))';[x x.^2 betainc(x,2,1)./x.^2-1]
ans =
1.000000000000000e+00 1.000000000000000e+00
0
5.000000000000000e-01 2.500000000000000e-01
2.220446049250313e-16
2.500000000000000e-01 6.250000000000000e-02
2.220446049250313e-16
1.250000000000000e-01 1.562500000000000e-02
2.220446049250313e-16
6.250000000000000e-02 3.906250000000000e-03
2.220446049250313e-16
3.125000000000000e-02 9.765625000000000e-04
-4.440892098500626e-16
1.562500000000000e-02 2.441406250000000e-04
8.881784197001252e-16
7.812500000000000e-03 6.103515625000000e-05
-1.110223024625157e-16
3.906250000000000e-03 1.525878906250000e-05
1.110223024625157e-15
1.953125000000000e-03 3.814697265625000e-06
-7.771561172376096e-16
9.765625000000000e-04 9.536743164062500e-07
-4.440892098500626e-16
4.882812500000000e-04 2.384185791015625e-07
4.440892098500626e-16
2.441406250000000e-04 5.960464477539063e-08
0
1.220703125000000e-04 1.490116119384766e-08
0
6.103515625000000e-05 3.725290298461914e-09
2.220446049250313e-15
3.051757812500000e-05 9.313225746154785e-10
-2.220446049250313e-16
1.525878906250000e-05 2.328306436538696e-10
1.776356839400250e-15
7.629394531250000e-06 5.820766091346741e-11
6.661338147750939e-16
3.814697265625000e-06 1.455191522836685e-11
4.440892098500626e-16
1.907348632812500e-06 3.637978807091713e-12
2.220446049250313e-16
9.536743164062500e-07 9.094947017729282e-13
0
>> x=2.^-18*(1+eps(1)*(0:9)');[x x.^2 betainc(x,2,1) betainc(x,2,1)./x.^2-1]
ans =
3.814697265625000e-06 1.455191522836685e-11 1.455191522836686e-11
4.440892098500626e-16
3.814697265625001e-06 1.455191522836686e-11 1.455191522836686e-11
0
3.814697265625002e-06 1.455191522836686e-11 1.455191522836686e-11
-4.440892098500626e-16
3.814697265625003e-06 1.455191522836687e-11 1.455191522836686e-11
-8.881784197001252e-16
3.814697265625003e-06 1.455191522836688e-11 1.455191522836686e-11
-1.332267629550188e-15
3.814697265625004e-06 1.455191522836688e-11 1.455191522836691e-11
1.776356839400250e-15
3.814697265625005e-06 1.455191522836689e-11 1.455191522836691e-11
1.332267629550188e-15
3.814697265625006e-06 1.455191522836690e-11 1.455191522836691e-11
8.881784197001252e-16
3.814697265625007e-06 1.455191522836690e-11 1.455191522836691e-11
4.440892098500626e-16
3.814697265625008e-06 1.455191522836691e-11 1.455191522836691e-11
0
using Octave 6.3.0:
>> x=2.^(-(0:20))';[x x.^2 betainc(x,2,1)./x.^2-1]
ans =
1.000000000000000e+00 1.000000000000000e+00 0
5.000000000000000e-01 2.500000000000000e-01 2.220446049250313e-16
2.500000000000000e-01 6.250000000000000e-02 0
1.250000000000000e-01 1.562500000000000e-02 4.440892098500626e-16
6.250000000000000e-02 3.906250000000000e-03 2.220446049250313e-16
3.125000000000000e-02 9.765625000000000e-04 0
1.562500000000000e-02 2.441406250000000e-04 6.661338147750939e-16
7.812500000000000e-03 6.103515625000000e-05 4.440892098500626e-16
3.906250000000000e-03 1.525878906250000e-05 4.440892098500626e-16
1.953125000000000e-03 3.814697265625000e-06 2.220446049250313e-16
9.765625000000000e-04 9.536743164062500e-07 0
4.882812500000000e-04 2.384185791015625e-07 -1.110223024625157e-16
2.441406250000000e-04 5.960464477539062e-08 1.554312234475219e-15
1.220703125000000e-04 1.490116119384766e-08 -2.331468351712829e-15
6.103515625000000e-05 3.725290298461914e-09 1.110223024625157e-15
3.051757812500000e-05 9.313225746154785e-10 -2.664535259100376e-15
1.525878906250000e-05 2.328306436538696e-10 6.661338147750939e-16
7.629394531250000e-06 5.820766091346741e-11 6.661338147750939e-16
3.814697265625000e-06 1.455191522836685e-11 -3.108624468950438e-15
1.907348632812500e-06 3.637978807091713e-12 2.220446049250313e-16
9.536743164062500e-07 9.094947017729282e-13 0
>> x=2.^-18*(1+eps(1)*(0:9)');[x x.^2 betainc(x,2,1) betainc(x,2,1)./x.^2-1]
ans =
3.814697265625000e-06 1.455191522836685e-11 1.455191522836681e-11
-3.108624468950438e-15
3.814697265625001e-06 1.455191522836686e-11 1.455191522836681e-11
-3.552713678800501e-15
3.814697265625002e-06 1.455191522836686e-11 1.455191522836681e-11
-3.996802888650564e-15
3.814697265625003e-06 1.455191522836687e-11 1.455191522836681e-11
-4.440892098500626e-15
3.814697265625003e-06 1.455191522836688e-11 1.455191522836681e-11
-4.884981308350689e-15
3.814697265625004e-06 1.455191522836688e-11 1.455191522836691e-11
1.776356839400250e-15
3.814697265625005e-06 1.455191522836689e-11 1.455191522836691e-11
1.332267629550188e-15
3.814697265625006e-06 1.455191522836690e-11 1.455191522836691e-11
8.881784197001252e-16
3.814697265625007e-06 1.455191522836690e-11 1.455191522836691e-11
4.440892098500626e-16
3.814697265625008e-06 1.455191522836691e-11 1.455191522836691e-11
0
I thought maybe the patch was already pushed but
Octave 6.3.0:
>> betainc (0.5, 1, Inf)
ans = NaN
after applying the patch i get:
>> x=2.^(-(0:20))';[x x.^2 betainc(x,2,1)./x.^2-1]
ans =
1.000000000000000e+00 1.000000000000000e+00 0
5.000000000000000e-01 2.500000000000000e-01 0
2.500000000000000e-01 6.250000000000000e-02 0
1.250000000000000e-01 1.562500000000000e-02 0
6.250000000000000e-02 3.906250000000000e-03 0
3.125000000000000e-02 9.765625000000000e-04 0
1.562500000000000e-02 2.441406250000000e-04 0
7.812500000000000e-03 6.103515625000000e-05 0
3.906250000000000e-03 1.525878906250000e-05 0
1.953125000000000e-03 3.814697265625000e-06 0
9.765625000000000e-04 9.536743164062500e-07 0
4.882812500000000e-04 2.384185791015625e-07 0
2.441406250000000e-04 5.960464477539062e-08 0
1.220703125000000e-04 1.490116119384766e-08 0
6.103515625000000e-05 3.725290298461914e-09 0
3.051757812500000e-05 9.313225746154785e-10 0
1.525878906250000e-05 2.328306436538696e-10 0
7.629394531250000e-06 5.820766091346741e-11 0
3.814697265625000e-06 1.455191522836685e-11 0
1.907348632812500e-06 3.637978807091713e-12 0
9.536743164062500e-07 9.094947017729282e-13 0
>> x=2.^-18*(1+eps(1)*(0:9)');[x x.^2 betainc(x,2,1) betainc(x,2,1)./x.^2-1]
ans =
3.814697265625000e-06 1.455191522836685e-11 1.455191522836685e-11
0
3.814697265625001e-06 1.455191522836686e-11 1.455191522836686e-11
0
3.814697265625002e-06 1.455191522836686e-11 1.455191522836686e-11
0
3.814697265625003e-06 1.455191522836687e-11 1.455191522836687e-11
0
3.814697265625003e-06 1.455191522836688e-11 1.455191522836688e-11
0
3.814697265625004e-06 1.455191522836688e-11 1.455191522836688e-11
0
3.814697265625005e-06 1.455191522836689e-11 1.455191522836689e-11
0
3.814697265625006e-06 1.455191522836690e-11 1.455191522836690e-11
0
3.814697265625007e-06 1.455191522836690e-11 1.455191522836690e-11
0
3.814697265625008e-06 1.455191522836691e-11 1.455191522836691e-11
0
>> betainc (0.5, 1, Inf)
ans = 1.000000000000000e+00
the trivial case fix is a good fix. but I'm not seeing the accuracy issue
_______________________________________________________
Reply to this item at:
<https://savannah.gnu.org/bugs/?60682>
_______________________________________________
Message sent via Savannah
https://savannah.gnu.org/
- [Octave-bug-tracker] [bug #60682] betainc is inaccurate,
Nicholas Jankowski <=