Compare commits

..

2 commits

Author SHA1 Message Date
Philip Dubé
cf775345f2
Merge 5dbfbaff9d into 2b68c56fe0 2024-04-01 23:05:16 +02:00
Mateusz Woś
2b68c56fe0
Adjust Ln method to prevent infinity iteration loops (#357)
* Adjust Ln method to prevent infinity iteration loops
* Add test case for infinity loop
2024-04-01 22:02:25 +02:00
2 changed files with 17 additions and 6 deletions

View file

@ -133,7 +133,6 @@ func NewFromBigInt(value *big.Int, exp int32) Decimal {
// d2 := NewFromBigRat(big.NewRat(4, 5), 1) // output: "0.8" // d2 := NewFromBigRat(big.NewRat(4, 5), 1) // output: "0.8"
// d3 := NewFromBigRat(big.NewRat(1000, 3), 3) // output: "333.333" // d3 := NewFromBigRat(big.NewRat(1000, 3), 3) // output: "333.333"
// d4 := NewFromBigRat(big.NewRat(2, 7), 4) // output: "0.2857" // d4 := NewFromBigRat(big.NewRat(2, 7), 4) // output: "0.2857"
//
func NewFromBigRat(value *big.Rat, precision int32) Decimal { func NewFromBigRat(value *big.Rat, precision int32) Decimal {
return Decimal{ return Decimal{
value: new(big.Int).Set(value.Num()), value: new(big.Int).Set(value.Num()),
@ -920,7 +919,10 @@ func (d Decimal) Ln(precision int32) (Decimal, error) {
// Halley's Iteration. // Halley's Iteration.
// Calculating n-th term of formula: a_(n+1) = a_n - 2 * (exp(a_n) - z) / (exp(a_n) + z), // Calculating n-th term of formula: a_(n+1) = a_n - 2 * (exp(a_n) - z) / (exp(a_n) + z),
// until the difference between current and next term is smaller than epsilon // until the difference between current and next term is smaller than epsilon
for { var prevStep Decimal
maxIters := calcPrecision*2 + 10
for i := int32(0); i < maxIters; i++ {
// exp(a_n) // exp(a_n)
comp3, _ = comp1.ExpTaylor(calcPrecision) comp3, _ = comp1.ExpTaylor(calcPrecision)
// exp(a_n) - z // exp(a_n) - z
@ -934,9 +936,17 @@ func (d Decimal) Ln(precision int32) (Decimal, error) {
// comp1 = a_(n+1) = a_n - 2 * (exp(a_n) - z) / (exp(a_n) + z) // comp1 = a_(n+1) = a_n - 2 * (exp(a_n) - z) / (exp(a_n) + z)
comp1 = comp1.Sub(comp3) comp1 = comp1.Sub(comp3)
if prevStep.Add(comp3).IsZero() {
// If iteration steps oscillate we should return early and prevent an infinity loop
// NOTE(mwoss): This should be quite a rare case, returning error is not necessary
break
}
if comp3.Abs().Cmp(epsilon) <= 0 { if comp3.Abs().Cmp(epsilon) <= 0 {
break break
} }
prevStep = comp3
} }
} }

View file

@ -2822,6 +2822,7 @@ func TestDecimal_Ln(t *testing.T) {
{"839101.0351094726488848490572028502", 50, "13.64008640145229044389152437468283605382056561604272"}, {"839101.0351094726488848490572028502", 50, "13.64008640145229044389152437468283605382056561604272"},
{"5023583755703750094849.03519358513093500275017501750602739169823", 25, "49.9684305274348922267409953"}, {"5023583755703750094849.03519358513093500275017501750602739169823", 25, "49.9684305274348922267409953"},
{"5023583755703750094849.03519358513093500275017501750602739169823", -1, "50.0"}, {"5023583755703750094849.03519358513093500275017501750602739169823", -1, "50.0"},
{"66.12", 18, "4.191471272952823429"},
} { } {
d, _ := NewFromString(testCase.Dec) d, _ := NewFromString(testCase.Dec)
expected, _ := NewFromString(testCase.Expected) expected, _ := NewFromString(testCase.Expected)