LEFT | RIGHT |
1 // Copyright 2010 The Go Authors. All rights reserved. | 1 // Copyright 2010 The Go Authors. All rights reserved. |
2 // Use of this source code is governed by a BSD-style | 2 // Use of this source code is governed by a BSD-style |
3 // license that can be found in the LICENSE file. | 3 // license that can be found in the LICENSE file. |
4 | 4 |
5 package big | 5 package big |
6 | 6 |
7 import ( | 7 import ( |
8 "bytes" | 8 "bytes" |
9 "encoding/gob" | 9 "encoding/gob" |
10 "fmt" | 10 "fmt" |
(...skipping 482 matching lines...) Expand 10 before | Expand all | Expand 10 after Loading... |
493 x.SetFrac(a, b) | 493 x.SetFrac(a, b) |
494 a = x.Num() | 494 a = x.Num() |
495 b = x.Denom() | 495 b = x.Denom() |
496 a.SetInt64(5) | 496 a.SetInt64(5) |
497 b.SetInt64(3) | 497 b.SetInt64(3) |
498 if x.Cmp(q53) != 0 { | 498 if x.Cmp(q53) != 0 { |
499 t.Errorf("3) got %s want %s", x, q53) | 499 t.Errorf("3) got %s want %s", x, q53) |
500 } | 500 } |
501 } | 501 } |
502 | 502 |
503 var inputs = []string{ | 503 // Test inputs to Rat.SetString. The optional prefix "slow:" skips |
| 504 // checks found to be slow for certain large rationals. |
| 505 var float64inputs = []string{ |
504 // | 506 // |
505 // Constants plundered from strconv/testfp.txt. | 507 // Constants plundered from strconv/testfp.txt. |
506 // | 508 // |
507 | 509 |
508 // Table 1: Stress Inputs for Conversion to 53-bit Binary, < 1/2 ULP | 510 // Table 1: Stress Inputs for Conversion to 53-bit Binary, < 1/2 ULP |
509 "5e+125", | 511 "5e+125", |
510 "69e+267", | 512 "69e+267", |
511 "999e-026", | 513 "999e-026", |
512 "7861e-034", | 514 "7861e-034", |
513 "75569e-254", | 515 "75569e-254", |
(...skipping 79 matching lines...) Expand 10 before | Expand all | Expand 10 after Loading... |
593 "100000000000000000000000", | 595 "100000000000000000000000", |
594 "1e-100", | 596 "1e-100", |
595 "123456700", | 597 "123456700", |
596 "99999999999999974834176", | 598 "99999999999999974834176", |
597 "100000000000000000000001", | 599 "100000000000000000000001", |
598 "100000000000000008388608", | 600 "100000000000000008388608", |
599 "100000000000000016777215", | 601 "100000000000000016777215", |
600 "100000000000000016777216", | 602 "100000000000000016777216", |
601 "-1", | 603 "-1", |
602 "-0.1", | 604 "-0.1", |
603 » "-0", | 605 » "-0", // NB: exception made for this input |
604 "1e-20", | 606 "1e-20", |
605 "625e-3", | 607 "625e-3", |
606 | 608 |
607 // largest float64 | 609 // largest float64 |
608 "1.7976931348623157e308", | 610 "1.7976931348623157e308", |
609 "-1.7976931348623157e308", | 611 "-1.7976931348623157e308", |
610 // next float64 - too large | 612 // next float64 - too large |
611 "1.7976931348623159e308", | 613 "1.7976931348623159e308", |
612 "-1.7976931348623159e308", | 614 "-1.7976931348623159e308", |
613 // the border is ...158079 | 615 // the border is ...158079 |
(...skipping 27 matching lines...) Expand all Loading... |
641 "1e-322", | 643 "1e-322", |
642 // smallest denormal | 644 // smallest denormal |
643 "5e-324", | 645 "5e-324", |
644 "4e-324", | 646 "4e-324", |
645 "3e-324", | 647 "3e-324", |
646 // too small | 648 // too small |
647 "2e-324", | 649 "2e-324", |
648 // way too small | 650 // way too small |
649 "1e-350", | 651 "1e-350", |
650 "slow:1e-400000", | 652 "slow:1e-400000", |
| 653 // way too small, negative |
| 654 "-1e-350", |
| 655 "slow:-1e-400000", |
651 | 656 |
652 // try to overflow exponent | 657 // try to overflow exponent |
653 // [Disabled: too slow and memory-hungry with rationals.] | 658 // [Disabled: too slow and memory-hungry with rationals.] |
654 // "1e-4294967296", | 659 // "1e-4294967296", |
655 // "1e+4294967296", | 660 // "1e+4294967296", |
656 // "1e-18446744073709551616", | 661 // "1e-18446744073709551616", |
657 // "1e+18446744073709551616", | 662 // "1e+18446744073709551616", |
658 | 663 |
659 // http://www.exploringbinary.com/java-hangs-when-converting-2-225073858
5072012e-308/ | 664 // http://www.exploringbinary.com/java-hangs-when-converting-2-225073858
5072012e-308/ |
660 "2.2250738585072012e-308", | 665 "2.2250738585072012e-308", |
(...skipping 21 matching lines...) Expand all Loading... |
682 // Smallest denormal, 2^(-1022-52) | 687 // Smallest denormal, 2^(-1022-52) |
683 "4.940656458412465441765687928682213723651e-324", | 688 "4.940656458412465441765687928682213723651e-324", |
684 // Half of smallest denormal, 2^(-1022-53) | 689 // Half of smallest denormal, 2^(-1022-53) |
685 "2.470328229206232720882843964341106861825e-324", | 690 "2.470328229206232720882843964341106861825e-324", |
686 // A little more than the exact half of smallest denormal | 691 // A little more than the exact half of smallest denormal |
687 // 2^-1075 + 2^-1100. (Rounds to 1p-1074.) | 692 // 2^-1075 + 2^-1100. (Rounds to 1p-1074.) |
688 "2.470328302827751011111470718709768633275e-324", | 693 "2.470328302827751011111470718709768633275e-324", |
689 // The exact halfway between smallest normal and largest denormal: | 694 // The exact halfway between smallest normal and largest denormal: |
690 // 2^-1022 - 2^-1075. (Rounds to 2^-1022.) | 695 // 2^-1022 - 2^-1075. (Rounds to 2^-1022.) |
691 "2.225073858507201136057409796709131975935e-308", | 696 "2.225073858507201136057409796709131975935e-308", |
| 697 |
| 698 "1152921504606846975", // 1<<60 - 1 |
| 699 "-1152921504606846975", // -(1<<60 - 1) |
| 700 "1152921504606846977", // 1<<60 + 1 |
| 701 "-1152921504606846977", // -(1<<60 + 1) |
| 702 |
| 703 "1/3", |
692 } | 704 } |
693 | 705 |
694 func TestFloat64SpecialCases(t *testing.T) { | 706 func TestFloat64SpecialCases(t *testing.T) { |
695 » for _, input := range inputs { | 707 » for _, input := range float64inputs { |
696 slow := strings.HasPrefix(input, "slow:") | 708 slow := strings.HasPrefix(input, "slow:") |
697 if slow { | 709 if slow { |
698 » » » input = input[5:] | 710 » » » input = input[len("slow:"):] |
699 } | 711 } |
700 | 712 |
701 r, ok := new(Rat).SetString(input) | 713 r, ok := new(Rat).SetString(input) |
702 if !ok { | 714 if !ok { |
703 » » » t.Error("can't parse rational", input) | 715 » » » t.Errorf("Rat.SetString(%q) failed", input) |
704 continue | 716 continue |
705 } | 717 } |
706 f, exact := r.Float64() | 718 f, exact := r.Float64() |
707 | 719 |
708 // 1. Check string -> Rat -> float64 conversions are | 720 // 1. Check string -> Rat -> float64 conversions are |
709 // consistent with strconv.ParseFloat. | 721 // consistent with strconv.ParseFloat. |
710 » » if e, _ := strconv.ParseFloat(input, 64); e != f { | 722 » » // Skip this check if the input uses "a/b" rational syntax. |
711 » » » t.Errorf("string -> Rat -> float64 yielded wrong value;
input=%s expected=%g actual=%g error=%g", input, e, f, (f-e)/e) | 723 » » if !strings.Contains(input, "/") { |
| 724 » » » e, _ := strconv.ParseFloat(input, 64) |
| 725 |
| 726 » » » // Careful: negative Rats too small for |
| 727 » » » // float64 become -0, but Rat obviously cannot |
| 728 » » » // preserve the sign from SetString("-0"). |
| 729 » » » switch { |
| 730 » » » case math.Float64bits(e) == math.Float64bits(f): |
| 731 » » » » // Ok: bitwise equal. |
| 732 » » » case f == 0 && r.Num().BitLen() == 0: |
| 733 » » » » // Ok: Rat(0) is equivalent to both +/- float64(
0). |
| 734 » » » default: |
| 735 » » » » t.Errorf("strconv.ParseFloat(%q) = %g (%b), want
%g (%b); delta=%g", input, e, e, f, f, f-e) |
| 736 » » » } |
712 } | 737 } |
713 | 738 |
714 if !isFinite(f) || slow { | 739 if !isFinite(f) || slow { |
715 continue | 740 continue |
716 } | 741 } |
717 | 742 |
718 // 2. Check f is best approximation to r. | 743 // 2. Check f is best approximation to r. |
719 if !checkIsBestApprox(t, f, r) { | 744 if !checkIsBestApprox(t, f, r) { |
720 » » » t.Errorf("Input was: %s", input) | 745 » » » // Append context information. |
| 746 » » » t.Errorf("(input was %q)", input) |
721 } | 747 } |
722 | 748 |
723 // 3. Check f->R->f roundtrip is non-lossy. | 749 // 3. Check f->R->f roundtrip is non-lossy. |
724 checkNonLossyRoundtrip(t, f) | 750 checkNonLossyRoundtrip(t, f) |
725 | 751 |
726 » » // 4. Check exactness. Slightly pointless since implementation
uses same algorithm. | 752 » » // 4. Check exactness using slow algorithm. |
727 if wasExact := new(Rat).SetFloat64(f).Cmp(r) == 0; wasExact != e
xact { | 753 if wasExact := new(Rat).SetFloat64(f).Cmp(r) == 0; wasExact != e
xact { |
728 » » » t.Errorf("Float64() returned exact=%b but actual exactne
ss was %b; input=%s", exact, wasExact, input) | 754 » » » t.Errorf("Rat.SetString(%q).Float64().exact = %b, want %
b", input, exact, wasExact) |
729 } | 755 } |
730 } | 756 } |
731 } | 757 } |
732 | 758 |
733 func TestFloat64Distribution(t *testing.T) { | 759 func TestFloat64Distribution(t *testing.T) { |
734 // Generate a distribution of (sign, mantissa, exp) values | 760 // Generate a distribution of (sign, mantissa, exp) values |
735 // broader than the float64 range, and check Rat.Float64() | 761 // broader than the float64 range, and check Rat.Float64() |
736 // always picks the closest float64 approximation. | 762 // always picks the closest float64 approximation. |
737 var add = []int64{ | 763 var add = []int64{ |
738 0, | 764 0, |
(...skipping 17 matching lines...) Expand all Loading... |
756 num, den := NewInt(b), NewInt(1) | 782 num, den := NewInt(b), NewInt(1) |
757 if exp > 0 { | 783 if exp > 0 { |
758 num.Lsh(num, uint(exp)) | 784 num.Lsh(num, uint(exp)) |
759 } else { | 785 } else { |
760 den.Lsh(den, uint(-exp)) | 786 den.Lsh(den, uint(-exp)) |
761 } | 787 } |
762 r := new(Rat).SetFrac(num, den) | 788 r := new(Rat).SetFrac(num, den) |
763 f, _ := r.Float64() | 789 f, _ := r.Float64() |
764 | 790 |
765 if !checkIsBestApprox(t, f, r) { | 791 if !checkIsBestApprox(t, f, r) { |
766 » » » » » » t.Errorf("Inputs were: %#x << %d
; f=%g; f~%g; r=%v", | 792 » » » » » » // Append context information. |
767 » » » » » » » b, exp, f, math.Ldexp(fl
oat64(b), exp), r) | 793 » » » » » » t.Errorf("(input was mantissa %#
x, exp %d; f=%g (%b); f~%g; r=%v)", |
| 794 » » » » » » » b, exp, f, f, math.Ldexp
(float64(b), exp), r) |
768 } | 795 } |
769 | 796 |
770 checkNonLossyRoundtrip(t, f) | 797 checkNonLossyRoundtrip(t, f) |
771 } | 798 } |
772 } | 799 } |
773 } | 800 } |
774 } | 801 } |
775 } | 802 } |
776 | 803 |
777 // TestFloat64NonFinite checks that SetFloat64 of a non-finite value | 804 // TestFloat64NonFinite checks that SetFloat64 of a non-finite value |
778 // returns nil. | 805 // returns nil. |
779 func TestSetFloat64NonFinite(t *testing.T) { | 806 func TestSetFloat64NonFinite(t *testing.T) { |
780 for _, f := range []float64{math.NaN(), math.Inf(+1), math.Inf(-1)} { | 807 for _, f := range []float64{math.NaN(), math.Inf(+1), math.Inf(-1)} { |
781 var r Rat | 808 var r Rat |
782 if r2 := r.SetFloat64(f); r2 != nil { | 809 if r2 := r.SetFloat64(f); r2 != nil { |
783 » » » t.Errorf("SetFloat64(%g) returned non-nil value: %v", f,
r2) | 810 » » » t.Errorf("SetFloat64(%g) was %v, want nil", f, r2) |
784 } | 811 } |
785 } | 812 } |
786 } | 813 } |
787 | 814 |
788 // checkNonLossyRoundtrip checks that a float->Rat->float roundtrip is | 815 // checkNonLossyRoundtrip checks that a float->Rat->float roundtrip is |
789 // non-lossy for finite f. | 816 // non-lossy for finite f. |
790 func checkNonLossyRoundtrip(t *testing.T, f float64) { | 817 func checkNonLossyRoundtrip(t *testing.T, f float64) { |
791 if !isFinite(f) { | 818 if !isFinite(f) { |
792 return | 819 return |
793 } | 820 } |
794 r := new(Rat).SetFloat64(f) | 821 r := new(Rat).SetFloat64(f) |
795 if r == nil { | 822 if r == nil { |
796 » » t.Error("SetFloat64 returned nil for a finite value", f) | 823 » » t.Errorf("Rat.SetFloat64(%g (%b)) == nil", f, f) |
797 return | 824 return |
798 } | 825 } |
799 f2, exact := r.Float64() | 826 f2, exact := r.Float64() |
800 » if f != f2 { | 827 » if f != f2 || !exact { |
801 » » t.Error("float64->Rat->float64 roundtrip was lossy", f, f2, (f2-
f)/f) | 828 » » t.Errorf("Rat.SetFloat64(%g).Float64() = %g (%b), %v, want %g (%
b), %v; delta=%b", |
802 » } | 829 » » » f, f2, f2, exact, f, f, true, f2-f) |
803 » if !exact { | 830 » } |
804 » » t.Error("float64->Rat->float64 roundtrip was inexact", f, f2) | 831 } |
805 » } | 832 |
| 833 // delta returns the absolute difference between r and f. |
| 834 func delta(r *Rat, f float64) *Rat { |
| 835 » d := new(Rat).Sub(r, new(Rat).SetFloat64(f)) |
| 836 » return d.Abs(d) |
806 } | 837 } |
807 | 838 |
808 // checkIsBestApprox checks that f is the best possible float64 | 839 // checkIsBestApprox checks that f is the best possible float64 |
809 // approximation of r. | 840 // approximation of r. |
810 // Returns true on success. | 841 // Returns true on success. |
811 func checkIsBestApprox(t *testing.T, f float64, r *Rat) bool { | 842 func checkIsBestApprox(t *testing.T, f float64, r *Rat) bool { |
812 if math.Abs(f) >= math.MaxFloat64 { | 843 if math.Abs(f) >= math.MaxFloat64 { |
813 // Cannot check +Inf, -Inf, nor the float next to them (MaxFloat
64). | 844 // Cannot check +Inf, -Inf, nor the float next to them (MaxFloat
64). |
814 // But we have tests for these special cases. | 845 // But we have tests for these special cases. |
815 return true | 846 return true |
816 } | 847 } |
817 | 848 |
818 // r must be strictly between f0 and f1, the floats bracketing f. | 849 // r must be strictly between f0 and f1, the floats bracketing f. |
819 f0 := math.Nextafter(f, math.Inf(-1)) | 850 f0 := math.Nextafter(f, math.Inf(-1)) |
820 f1 := math.Nextafter(f, math.Inf(+1)) | 851 f1 := math.Nextafter(f, math.Inf(+1)) |
821 » rf := new(Rat).SetFloat64(f) // the approximation as a Rat | 852 |
822 » rf0 := new(Rat).SetFloat64(f0) | 853 » // For f to be correct, r must be closer to f than to f0 or f1. |
823 » rf1 := new(Rat).SetFloat64(f1) | 854 » df := delta(r, f) |
824 » if r.Cmp(rf0) <= 0 || r.Cmp(rf1) >= 0 { | 855 » df0 := delta(r, f0) |
825 » » t.Errorf("%b (approx. of %v) is not in interval [%b, %b]", f, r,
f0, f1) | 856 » df1 := delta(r, f1) |
| 857 » if df.Cmp(df0) > 0 { |
| 858 » » t.Errorf("Rat(%v).Float64() = %g (%b), but previous float64 %g (
%b) is closer", r, f, f, f0, f0) |
826 return false | 859 return false |
827 } | 860 } |
828 | 861 » if df.Cmp(df1) > 0 { |
829 » // For f to be correct, r must be closer to f than to f0 or f1. | 862 » » t.Errorf("Rat(%v).Float64() = %g (%b), but next float64 %g (%b)
is closer", r, f, f, f1, f1) |
830 » df := new(Rat).Sub(r, rf).Abs(new(Rat)) | |
831 » df0 := new(Rat).Sub(r, rf0) // known > 0 | |
832 » df1 := new(Rat).Sub(rf1, r) // known > 0 | |
833 » if df.Cmp(df0) > 0 { | |
834 » » t.Errorf("%b (approx. of %v) is not as close as previous float %
b", f, r, f0) | |
835 return false | 863 return false |
836 } | 864 } |
837 » if df.Cmp(df1) > 0 { | 865 » if df.Cmp(df0) == 0 && !isEven(f) { |
838 » » t.Errorf("%b (approx. of %v) is not as close as next float %b",
f, r, f1) | 866 » » t.Errorf("Rat(%v).Float64() = %g (%b); halfway should have round
ed to %g (%b) instead", r, f, f, f0, f0) |
839 return false | 867 return false |
840 } | 868 } |
841 if df.Cmp(df0) == 0 && !isEven(f) { | |
842 t.Errorf("%b (approx. of %v) should round halfway to previous (e
ven) float %b ", f, r, f0) | |
843 return false | |
844 } | |
845 if df.Cmp(df1) == 0 && !isEven(f) { | 869 if df.Cmp(df1) == 0 && !isEven(f) { |
846 » » t.Errorf("%b (approx. of %v) should round halfway to next (even)
float %b", f, r, f1) | 870 » » t.Errorf("Rat(%v).Float64() = %g (%b); halfway should have round
ed to %g (%b) instead", r, f, f, f1, f1) |
847 return false | 871 return false |
848 } | 872 } |
849 return true | 873 return true |
850 } | 874 } |
851 | 875 |
852 func isEven(f float64) bool { return math.Float64bits(f)&1 == 0 } | 876 func isEven(f float64) bool { return math.Float64bits(f)&1 == 0 } |
853 | 877 |
854 func TestIsFinite(t *testing.T) { | 878 func TestIsFinite(t *testing.T) { |
855 finites := []float64{ | 879 finites := []float64{ |
856 » » 1 / 3, | 880 » » 1.0 / 3, |
857 4891559871276714924261e+222, | 881 4891559871276714924261e+222, |
858 math.MaxFloat64, | 882 math.MaxFloat64, |
859 math.SmallestNonzeroFloat64, | 883 math.SmallestNonzeroFloat64, |
860 -math.MaxFloat64, | 884 -math.MaxFloat64, |
861 -math.SmallestNonzeroFloat64, | 885 -math.SmallestNonzeroFloat64, |
862 } | 886 } |
863 for _, f := range finites { | 887 for _, f := range finites { |
864 if !isFinite(f) { | 888 if !isFinite(f) { |
865 » » » t.Error("!IsFinite(", f, ")") | 889 » » » t.Errorf("!IsFinite(%g (%b))", f, f) |
866 } | 890 } |
867 } | 891 } |
868 nonfinites := []float64{ | 892 nonfinites := []float64{ |
869 math.NaN(), | 893 math.NaN(), |
870 math.Inf(-1), | 894 math.Inf(-1), |
871 math.Inf(+1), | 895 math.Inf(+1), |
872 } | 896 } |
873 for _, f := range nonfinites { | 897 for _, f := range nonfinites { |
874 if isFinite(f) { | 898 if isFinite(f) { |
875 » » » t.Error("IsFinite(", f, ")") | 899 » » » t.Errorf("IsFinite(%g, (%b))", f, f) |
876 » » } | 900 » » } |
877 » } | 901 » } |
878 } | 902 } |
LEFT | RIGHT |