Skip to content

Commit

Permalink
feat(osmomath): log2 approximation (#2788)
Browse files Browse the repository at this point in the history
* feat: osmomath log2 approximation

* lint

* fix comment

* comment

* Update osmomath/decimal.go

* improve accuracy with narrower range

* comment

* implement and test precise log for x >= 1

* tests for negative values

* make non-mutative and test

* bench

* changelog

* remove redundant assignments

* improve comments
  • Loading branch information
p0mvn authored Oct 24, 2022
1 parent 7119419 commit 8590b80
Show file tree
Hide file tree
Showing 5 changed files with 205 additions and 31 deletions.
1 change: 1 addition & 0 deletions CHANGELOG.md
Original file line number Diff line number Diff line change
Expand Up @@ -44,6 +44,7 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0

### Features

* [#2788](https://github.com/osmosis-labs/osmosis/pull/2788) Add logarithm base 2 implementation.
* [#2739](https://github.com/osmosis-labs/osmosis/pull/2739) Add pool type query

### Bug fixes
Expand Down
63 changes: 62 additions & 1 deletion osmomath/decimal.go
Original file line number Diff line number Diff line change
Expand Up @@ -30,6 +30,9 @@ const (

// max number of iterations in ApproxRoot function
maxApproxRootIterations = 100

// max number of iterations in Log2 function
maxLog2Iterations = 300
)

var (
Expand All @@ -39,6 +42,10 @@ var (
zeroInt = big.NewInt(0)
oneInt = big.NewInt(1)
tenInt = big.NewInt(10)

// initialized in init() since requires
// precision to be defined.
twoBigDec BigDec
)

// Decimal errors
Expand All @@ -54,6 +61,8 @@ func init() {
for i := 0; i <= Precision; i++ {
precisionMultipliers[i] = calcPrecisionMultiplier(int64(i))
}

twoBigDec = NewBigDec(2)
}

func precisionInt() *big.Int {
Expand Down Expand Up @@ -211,7 +220,8 @@ func (d BigDec) GTE(d2 BigDec) bool { return (d.i).Cmp(d2.i) >= 0 } /
func (d BigDec) LT(d2 BigDec) bool { return (d.i).Cmp(d2.i) < 0 } // less than
func (d BigDec) LTE(d2 BigDec) bool { return (d.i).Cmp(d2.i) <= 0 } // less than or equal
func (d BigDec) Neg() BigDec { return BigDec{new(big.Int).Neg(d.i)} } // reverse the decimal sign
func (d BigDec) Abs() BigDec { return BigDec{new(big.Int).Abs(d.i)} } // absolute value
// nolint: stylecheck
func (d BigDec) Abs() BigDec { return BigDec{new(big.Int).Abs(d.i)} } // absolute value

// BigInt returns a copy of the underlying big.Int.
func (d BigDec) BigInt() *big.Int {
Expand Down Expand Up @@ -856,3 +866,54 @@ func DecApproxEq(t *testing.T, d1 BigDec, d2 BigDec, tol BigDec) (*testing.T, bo
diff := d1.Sub(d2).Abs()
return t, diff.LTE(tol), "expected |d1 - d2| <:\t%v\ngot |d1 - d2| = \t\t%v", tol.String(), diff.String()
}

// LogBase2 returns log_2 {x}.
// Rounds down by truncations during division and right shifting.
// Accurate up to 32 precision digits.
// Implementation is based on:
// https://stm32duinoforum.com/forum/dsp/BinaryLogarithm.pdf
func (x BigDec) LogBase2() BigDec {
// create a new decimal to avoid mutating
// the receiver's int buffer.
xCopy := ZeroDec()
xCopy.i = new(big.Int).Set(x.i)
if xCopy.LTE(ZeroDec()) {
panic(fmt.Sprintf("log is not defined at <= 0, given (%s)", xCopy))
}

// Normalize x to be 1 <= x < 2.

// y is the exponent that results in a whole multiple of 2.
y := ZeroDec()

// repeat until: x >= 1.
for xCopy.LT(OneDec()) {
xCopy.i.Lsh(xCopy.i, 1)
y = y.Sub(OneDec())
}

// repeat until: x < 2.
for xCopy.GTE(twoBigDec) {
xCopy.i.Rsh(xCopy.i, 1)
y = y.Add(OneDec())
}

b := OneDec().Quo(twoBigDec)

// N.B. At this point x is a positive real number representing
// mantissa of the log. We estimate it using the following
// algorithm:
// https://stm32duinoforum.com/forum/dsp/BinaryLogarithm.pdf
// This has shown precision of 32 digits relative
// to Wolfram Alpha in tests.
for i := 0; i < maxLog2Iterations; i++ {
xCopy = xCopy.Mul(xCopy)
if xCopy.GTE(twoBigDec) {
xCopy.i.Rsh(xCopy.i, 1)
y = y.Add(b)
}
b.i.Rsh(b.i, 1)
}

return y
}
147 changes: 118 additions & 29 deletions osmomath/decimal_test.go
Original file line number Diff line number Diff line change
Expand Up @@ -11,6 +11,8 @@ import (
"github.com/stretchr/testify/require"
"github.com/stretchr/testify/suite"
"gopkg.in/yaml.v2"

"github.com/osmosis-labs/osmosis/v12/app/apptesting/osmoassert"
)

type decimalTestSuite struct {
Expand Down Expand Up @@ -152,8 +154,8 @@ func (s *decimalTestSuite) TestDecFloat64() {

func (s *decimalTestSuite) TestSdkDec() {
tests := []struct {
d BigDec
want sdk.Dec
d BigDec
want sdk.Dec
expPanic bool
}{
{NewBigDec(0), sdk.MustNewDecFromStr("0.000000000000000000"), false},
Expand All @@ -177,8 +179,8 @@ func (s *decimalTestSuite) TestSdkDec() {

func (s *decimalTestSuite) TestBigDecFromSdkDec() {
tests := []struct {
d sdk.Dec
want BigDec
d sdk.Dec
want BigDec
expPanic bool
}{
{sdk.MustNewDecFromStr("0.000000000000000000"), NewBigDec(0), false},
Expand All @@ -202,8 +204,8 @@ func (s *decimalTestSuite) TestBigDecFromSdkDec() {

func (s *decimalTestSuite) TestBigDecFromSdkDecSlice() {
tests := []struct {
d []sdk.Dec
want []BigDec
d []sdk.Dec
want []BigDec
expPanic bool
}{
{[]sdk.Dec{sdk.MustNewDecFromStr("0.000000000000000000")}, []BigDec{NewBigDec(0)}, false},
Expand Down Expand Up @@ -441,14 +443,14 @@ func (s *decimalTestSuite) TestDecCeil() {
input BigDec
expected BigDec
}{
{MustNewDecFromStr("0.001"), NewBigDec(1)}, // 0.001 => 1.0
{MustNewDecFromStr("-0.001"), ZeroDec()}, // -0.001 => 0.0
{ZeroDec(), ZeroDec()}, // 0.0 => 0.0
{MustNewDecFromStr("0.9"), NewBigDec(1)}, // 0.9 => 1.0
{MustNewDecFromStr("0.001"), NewBigDec(1)}, // 0.001 => 1.0
{MustNewDecFromStr("-0.001"), ZeroDec()}, // -0.001 => 0.0
{ZeroDec(), ZeroDec()}, // 0.0 => 0.0
{MustNewDecFromStr("0.9"), NewBigDec(1)}, // 0.9 => 1.0
{MustNewDecFromStr("4.001"), NewBigDec(5)}, // 4.001 => 5.0
{MustNewDecFromStr("-4.001"), NewBigDec(-4)}, // -4.001 => -4.0
{MustNewDecFromStr("4.7"), NewBigDec(5)}, // 4.7 => 5.0
{MustNewDecFromStr("-4.7"), NewBigDec(-4)}, // -4.7 => -4.0
{MustNewDecFromStr("4.7"), NewBigDec(5)}, // 4.7 => 5.0
{MustNewDecFromStr("-4.7"), NewBigDec(-4)}, // -4.7 => -4.0
}

for i, tc := range testCases {
Expand All @@ -463,11 +465,11 @@ func (s *decimalTestSuite) TestPower() {
power uint64
expected BigDec
}{
{OneDec(), 10, OneDec()}, // 1.0 ^ (10) => 1.0
{NewDecWithPrec(5, 1), 2, NewDecWithPrec(25, 2)}, // 0.5 ^ 2 => 0.25
{NewDecWithPrec(2, 1), 2, NewDecWithPrec(4, 2)}, // 0.2 ^ 2 => 0.04
{NewDecFromInt(NewInt(3)), 3, NewDecFromInt(NewInt(27))}, // 3 ^ 3 => 27
{NewDecFromInt(NewInt(-3)), 4, NewDecFromInt(NewInt(81))}, // -3 ^ 4 = 81
{OneDec(), 10, OneDec()}, // 1.0 ^ (10) => 1.0
{NewDecWithPrec(5, 1), 2, NewDecWithPrec(25, 2)}, // 0.5 ^ 2 => 0.25
{NewDecWithPrec(2, 1), 2, NewDecWithPrec(4, 2)}, // 0.2 ^ 2 => 0.04
{NewDecFromInt(NewInt(3)), 3, NewDecFromInt(NewInt(27))}, // 3 ^ 3 => 27
{NewDecFromInt(NewInt(-3)), 4, NewDecFromInt(NewInt(81))}, // -3 ^ 4 = 81
{MustNewDecFromStr("1.414213562373095048801688724209698079"), 2, NewDecFromInt(NewInt(2))}, // 1.414213562373095048801688724209698079 ^ 2 = 2
}

Expand All @@ -483,14 +485,14 @@ func (s *decimalTestSuite) TestApproxRoot() {
root uint64
expected BigDec
}{
{OneDec(), 10, OneDec()}, // 1.0 ^ (0.1) => 1.0
{NewDecWithPrec(25, 2), 2, NewDecWithPrec(5, 1)}, // 0.25 ^ (0.5) => 0.5
{NewDecWithPrec(4, 2), 2, NewDecWithPrec(2, 1)}, // 0.04 ^ (0.5) => 0.2
{NewDecFromInt(NewInt(27)), 3, NewDecFromInt(NewInt(3))}, // 27 ^ (1/3) => 3
{NewDecFromInt(NewInt(-81)), 4, NewDecFromInt(NewInt(-3))}, // -81 ^ (0.25) => -3
{NewDecFromInt(NewInt(2)), 2, MustNewDecFromStr("1.414213562373095048801688724209698079")}, // 2 ^ (0.5) => 1.414213562373095048801688724209698079
{OneDec(), 10, OneDec()}, // 1.0 ^ (0.1) => 1.0
{NewDecWithPrec(25, 2), 2, NewDecWithPrec(5, 1)}, // 0.25 ^ (0.5) => 0.5
{NewDecWithPrec(4, 2), 2, NewDecWithPrec(2, 1)}, // 0.04 ^ (0.5) => 0.2
{NewDecFromInt(NewInt(27)), 3, NewDecFromInt(NewInt(3))}, // 27 ^ (1/3) => 3
{NewDecFromInt(NewInt(-81)), 4, NewDecFromInt(NewInt(-3))}, // -81 ^ (0.25) => -3
{NewDecFromInt(NewInt(2)), 2, MustNewDecFromStr("1.414213562373095048801688724209698079")}, // 2 ^ (0.5) => 1.414213562373095048801688724209698079
{NewDecWithPrec(1005, 3), 31536000, MustNewDecFromStr("1.000000000158153903837946258002096839")}, // 1.005 ^ (1/31536000) ≈ 1.000000000158153903837946258002096839
{SmallestDec(), 2, NewDecWithPrec(1, 18)}, // 1e-36 ^ (0.5) => 1e-18
{SmallestDec(), 2, NewDecWithPrec(1, 18)}, // 1e-36 ^ (0.5) => 1e-18
{SmallestDec(), 3, MustNewDecFromStr("0.000000000001000000000000000002431786")}, // 1e-36 ^ (1/3) => 1e-12
{NewDecWithPrec(1, 8), 3, MustNewDecFromStr("0.002154434690031883721759293566519280")}, // 1e-8 ^ (1/3) ≈ 0.002154434690031883721759293566519
}
Expand All @@ -511,11 +513,11 @@ func (s *decimalTestSuite) TestApproxSqrt() {
input BigDec
expected BigDec
}{
{OneDec(), OneDec()}, // 1.0 => 1.0
{NewDecWithPrec(25, 2), NewDecWithPrec(5, 1)}, // 0.25 => 0.5
{NewDecWithPrec(4, 2), NewDecWithPrec(2, 1)}, // 0.09 => 0.3
{NewDecFromInt(NewInt(9)), NewDecFromInt(NewInt(3))}, // 9 => 3
{NewDecFromInt(NewInt(-9)), NewDecFromInt(NewInt(-3))}, // -9 => -3
{OneDec(), OneDec()}, // 1.0 => 1.0
{NewDecWithPrec(25, 2), NewDecWithPrec(5, 1)}, // 0.25 => 0.5
{NewDecWithPrec(4, 2), NewDecWithPrec(2, 1)}, // 0.09 => 0.3
{NewDecFromInt(NewInt(9)), NewDecFromInt(NewInt(3))}, // 9 => 3
{NewDecFromInt(NewInt(-9)), NewDecFromInt(NewInt(-3))}, // -9 => -3
{NewDecFromInt(NewInt(2)), MustNewDecFromStr("1.414213562373095048801688724209698079")}, // 2 => 1.414213562373095048801688724209698079
}

Expand Down Expand Up @@ -650,3 +652,90 @@ func BenchmarkMarshalTo(b *testing.B) {
}
}
}

func (s *decimalTestSuite) TestLog2() {
var expectedErrTolerance = MustNewDecFromStr("0.000000000000000000000000000000000100")

tests := map[string]struct {
initialValue BigDec
expected BigDec

expectedPanic bool
}{
"log_2{-1}; invalid; panic": {
initialValue: OneDec().Neg(),
expectedPanic: true,
},
"log_2{0}; invalid; panic": {
initialValue: ZeroDec(),
expectedPanic: true,
},
"log_2{0.001} = -9.965784284662087043610958288468170528": {
initialValue: MustNewDecFromStr("0.001"),
// From: https://www.wolframalpha.com/input?i=log+base+2+of+0.999912345+with+33+digits
expected: MustNewDecFromStr("-9.965784284662087043610958288468170528"),
},
"log_2{0.56171821941421412902170941} = -0.832081497183140708984033250637831402": {
initialValue: MustNewDecFromStr("0.56171821941421412902170941"),
// From: https://www.wolframalpha.com/input?i=log+base+2+of+0.56171821941421412902170941+with+36+digits
expected: MustNewDecFromStr("-0.832081497183140708984033250637831402"),
},
"log_2{0.999912345} = -0.000126464976533858080645902722235833": {
initialValue: MustNewDecFromStr("0.999912345"),
// From: https://www.wolframalpha.com/input?i=log+base+2+of+0.999912345+with+37+digits
expected: MustNewDecFromStr("-0.000126464976533858080645902722235833"),
},
"log_2{1} = 0": {
initialValue: NewBigDec(1),
expected: NewBigDec(0),
},
"log_2{2} = 1": {
initialValue: NewBigDec(2),
expected: NewBigDec(1),
},
"log_2{7} = 2.807354922057604107441969317231830809": {
initialValue: NewBigDec(7),
// From: https://www.wolframalpha.com/input?i=log+base+2+of+7+37+digits
expected: MustNewDecFromStr("2.807354922057604107441969317231830809"),
},
"log_2{512} = 9": {
initialValue: NewBigDec(512),
expected: NewBigDec(9),
},
"log_2{580} = 9.179909090014934468590092754117374938": {
initialValue: NewBigDec(580),
// From: https://www.wolframalpha.com/input?i=log+base+2+of+600+37+digits
expected: MustNewDecFromStr("9.179909090014934468590092754117374938"),
},
"log_2{1024} = 10": {
initialValue: NewBigDec(1024),
expected: NewBigDec(10),
},
"log_2{1024.987654321} = 10.001390817654141324352719749259888355": {
initialValue: NewDecWithPrec(1024987654321, 9),
// From: https://www.wolframalpha.com/input?i=log+base+2+of+1024.987654321+38+digits
expected: MustNewDecFromStr("10.001390817654141324352719749259888355"),
},
"log_2{912648174127941279170121098210.92821920190204131121} = 99.525973560175362367047484597337715868": {
initialValue: MustNewDecFromStr("912648174127941279170121098210.92821920190204131121"),
// From: https://www.wolframalpha.com/input?i=log+base+2+of+912648174127941279170121098210.92821920190204131121+38+digits
expected: MustNewDecFromStr("99.525973560175362367047484597337715868"),
},
}

for name, tc := range tests {
s.Run(name, func() {
osmoassert.ConditionalPanic(s.T(), tc.expectedPanic, func() {
// Create a copy to test that the original was not modified.
// That is, that LogbBase2() is non-mutative.
initialCopy := ZeroDec()
initialCopy.i.Set(tc.initialValue.i)

// system under test.
res := tc.initialValue.LogBase2()
require.True(DecApproxEq(s.T(), tc.expected, res, expectedErrTolerance))
require.Equal(s.T(), initialCopy, tc.initialValue)
})
})
}
}
24 changes: 24 additions & 0 deletions osmomath/log2_bench_test.go
Original file line number Diff line number Diff line change
@@ -0,0 +1,24 @@
package osmomath

import (
"math/rand"
"testing"
)

func BenchmarkLog2(b *testing.B) {
tests := []BigDec{
MustNewDecFromStr("1.2"),
MustNewDecFromStr("1.234"),
MustNewDecFromStr("1024"),
NewBigDec(2048 * 2048 * 2048 * 2048 * 2048),
MustNewDecFromStr("999999999999999999999999999999999999999999999999999999.9122181273612911"),
MustNewDecFromStr("0.563289239121902491248219047129047129"),
}

for i := 0; i < b.N; i++ {
b.StopTimer()
test := tests[rand.Int63n(int64(len(tests)))]
b.StartTimer()
_ = test.LogBase2()
}
}
1 change: 0 additions & 1 deletion osmomath/rounding_direction_test.go
Original file line number Diff line number Diff line change
Expand Up @@ -42,7 +42,6 @@ func TestDivIntByU64ToBigDec(t *testing.T) {
addTCForAllRoundingModes("odd divided by 2", sdk.NewInt(5), 2, NewDecWithPrec(25, 1))

for name, tt := range tests {
fmt.Println("start")
t.Run(name, func(t *testing.T) {
got, err := DivIntByU64ToBigDec(tt.i, tt.u, tt.round)
require.Equal(t, tt.want, got)
Expand Down

0 comments on commit 8590b80

Please sign in to comment.