Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

feat(osmomath): Exp2 function #3708

Merged
merged 20 commits into from
Dec 16, 2022
Merged
Show file tree
Hide file tree
Changes from 12 commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
1 change: 1 addition & 0 deletions CHANGELOG.md
Original file line number Diff line number Diff line change
Expand Up @@ -49,6 +49,7 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
* [#3677](https://github.com/osmosis-labs/osmosis/pull/3677) Add methods for cloning and mutative multiplication on osmomath.BigDec.
* [#3676](https://github.com/osmosis-labs/osmosis/pull/3676) implement `PowerInteger` function on `osmomath.BigDec`
* [#3678](https://github.com/osmosis-labs/osmosis/pull/3678) implement mutative `PowerIntegerMut` function on `osmomath.BigDec`.
* [#3708](https://github.com/osmosis-labs/osmosis/pull/3708) `Exp2` function to compute 2^decimal.
* [#3693](https://github.com/osmosis-labs/osmosis/pull/3693) Add `EstimateSwapExactAmountOut` query to stargate whitelist

### Bug fixes
Expand Down
92 changes: 92 additions & 0 deletions osmomath/exp2.go
Original file line number Diff line number Diff line change
@@ -0,0 +1,92 @@
package osmomath

import "fmt"

var (
// Truncated at precision end.
// See scripts/approximations/main.py exponent_approximation_choice function for details.
numeratorCoefficients13Param = []BigDec{
MustNewDecFromStr("1.000000000000000000000044212244679434"),
MustNewDecFromStr("0.352032455817400196452603772766844426"),
MustNewDecFromStr("0.056507868883666405413116800969512484"),
MustNewDecFromStr("0.005343900728213034434757419480319916"),
MustNewDecFromStr("0.000317708814342353603087543715930732"),
MustNewDecFromStr("0.000011429747507407623028722262874632"),
MustNewDecFromStr("0.000000198381965651614980168744540366"),
}

// Rounded up at precision end.
// See scripts/approximations/main.py exponent_approximation_choice function for details.
denominatorCoefficients13Param = []BigDec{
OneDec(),
MustNewDecFromStr("0.341114724742545112949699755780593311").Neg(),
MustNewDecFromStr("0.052724071627342653404436933178482287"),
MustNewDecFromStr("0.004760950735524957576233524801866342").Neg(),
MustNewDecFromStr("0.000267168475410566529819971616894193"),
MustNewDecFromStr("0.000008923715368802211181557353097439").Neg(),
MustNewDecFromStr("0.000000140277233177373698516010555916"),
}
p0mvn marked this conversation as resolved.
Show resolved Hide resolved

// maxSupportedExponent = 2^10. The value is chosen by benchmarking
// when the underlying internal functions overflow.
// If needed in the future, Exp2 can be reimplemented to allow for greater exponents.
maxSupportedExponent = MustNewDecFromStr("2").PowerInteger(9)
)

// Exp2 takes 2 to the power of a given non-negative decimal exponent
// and returns the result.
// The computation is performed by using th following property:
// 2^decimal_exp = 2^{integer_exp + fractional_exp} = 2^integer_exp * 2^fractional_exp
// The max supported exponent is defined by the global maxSupportedExponent.
// If a greater exponent is given, the function panics.
// Panics if the exponent is negative.
func Exp2(exponent BigDec) BigDec {
if exponent.Abs().GT(maxSupportedExponent) {
panic(fmt.Sprintf("integer exponent %s is too large, max (%s)", exponent, maxSupportedExponent))
}
if exponent.IsNegative() {
panic(fmt.Sprintf("negative exponent %s is not supported", exponent))
}

integerExponent := exponent.TruncateDec()

fractionalExponent := exponent.Sub(integerExponent)
fractionalResult := exp2ChebyshevRationalApprox(fractionalExponent)

fractionalResult.i = fractionalResult.i.Lsh(fractionalResult.i, uint(integerExponent.TruncateInt().Uint64()))

return fractionalResult
}

// exp2ChebyshevRationalApprox takes 2 to the power of a given decimal exponent.
// The result is approximated by a 13 parameter Chebyshev rational approximation.
// f(x) = h(x) / p(x) (7, 7) terms. We set the first term of p(x) to 1.
// As a result, this ends up being 7 + 6 = 13 parameters.
// The numerator coefficients are truncated at precision end. The denominator
// coefficients are rounded up at precision end.
// See scripts/approximations/README.md for details of the scripts used
// to compute the coefficients.
// CONTRACT: exponent must be in the range [0, 1], panics if not.
func exp2ChebyshevRationalApprox(x BigDec) BigDec {
if x.LT(ZeroDec()) || x.GT(OneDec()) {
panic(fmt.Sprintf("exponent must be in the range [0, 1], got %s", x))
}
if x.IsZero() {
return OneDec()
}
if x.Equal(OneDec()) {
return twoBigDec
}

h_x := numeratorCoefficients13Param[0].Clone()
p_x := denominatorCoefficients13Param[0].Clone()
x_exp_i := OneDec()
for i := 1; i < len(numeratorCoefficients13Param); i++ {
x_exp_i.MulMut(x)

h_x = h_x.Add(numeratorCoefficients13Param[i].Mul(x_exp_i))
p_x = p_x.Add(denominatorCoefficients13Param[i].Mul(x_exp_i))
}

return h_x.Quo(p_x)
p0mvn marked this conversation as resolved.
Show resolved Hide resolved
}
296 changes: 296 additions & 0 deletions osmomath/exp2_test.go
Original file line number Diff line number Diff line change
@@ -0,0 +1,296 @@
package osmomath_test

import (
"testing"

sdk "github.com/cosmos/cosmos-sdk/types"
"github.com/stretchr/testify/require"

"github.com/osmosis-labs/osmosis/v13/app/apptesting/osmoassert"
"github.com/osmosis-labs/osmosis/v13/osmomath"
"github.com/osmosis-labs/osmosis/v13/osmoutils"
)

var (
// minDecTolerance minimum tolerance for sdk.Dec, given its precision of 18.
minDecTolerance = sdk.MustNewDecFromStr("0.000000000000000001")
)

func TestExp2ChebyshevRationalApprox(t *testing.T) {
smallValue := osmomath.MustNewDecFromStr("0.00001")
smallerValue := osmomath.MustNewDecFromStr("0.00000000000000000001")
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

What is this smaller than & what is relation to minDecTolerance (sorry not sure if theres some important relation)

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I used these to implement test cases closer to the boundaries.

The reason is that the errors tend to be large at the boundaries for non-Chebyshev approximations. This is just to make sure that, as expected, Chebyshev does not allow for larger errors close to boundaries. I will add a comment relaying this


tests := map[string]struct {
exponent osmomath.BigDec
expectedResult osmomath.BigDec
errTolerance osmoutils.ErrTolerance
expectPanic bool
}{
"exp2(0.5)": {
exponent: osmomath.MustNewDecFromStr("0.5"),
// https://www.wolframalpha.com/input?i=2%5E0.5+37+digits
expectedResult: osmomath.MustNewDecFromStr("1.414213562373095048801688724209698079"),

errTolerance: osmoutils.ErrTolerance{
AdditiveTolerance: minDecTolerance,
MultiplicativeTolerance: minDecTolerance,
RoundingDir: osmomath.RoundDown,
},
},
"exp2(0)": {
exponent: osmomath.ZeroDec(),
expectedResult: osmomath.OneDec(),

errTolerance: osmoutils.ErrTolerance{
AdditiveTolerance: sdk.ZeroDec(),
MultiplicativeTolerance: sdk.ZeroDec(),
RoundingDir: osmomath.RoundDown,
},
},
"exp2(1)": {
exponent: osmomath.OneDec(),
expectedResult: osmomath.MustNewDecFromStr("2"),

errTolerance: osmoutils.ErrTolerance{
AdditiveTolerance: sdk.ZeroDec(),
MultiplicativeTolerance: sdk.ZeroDec(),
RoundingDir: osmomath.RoundDown,
},
},
"exp2(0.00001)": {
exponent: smallValue,
// https://www.wolframalpha.com/input?i=2%5E0.00001+37+digits
expectedResult: osmomath.MustNewDecFromStr("1.000006931495828305653209089800561681"),

errTolerance: osmoutils.ErrTolerance{
AdditiveTolerance: minDecTolerance,
MultiplicativeTolerance: minDecTolerance,
// TODO: confirm if rounding behavior is acceptable.
// Note, that a Python estimate passes RoundDown but not Wolfram.
RoundingDir: osmomath.RoundUnconstrained,
},
},
"exp2(0.99999)": {
exponent: osmomath.OneDec().Sub(smallValue),
// https://www.wolframalpha.com/input?i=2%5E0.99999+37+digits
expectedResult: osmomath.MustNewDecFromStr("1.999986137104433991477606830496602898"),

errTolerance: osmoutils.ErrTolerance{
AdditiveTolerance: sdk.MustNewDecFromStr("0.00000000000000007"),
MultiplicativeTolerance: minDecTolerance.Mul(sdk.NewDec(100)),
RoundingDir: osmomath.RoundDown,
},
},
"exp2(0.99999...)": {
exponent: osmomath.OneDec().Sub(smallerValue),
// https://www.wolframalpha.com/input?i=2%5E%281+-+0.00000000000000000001%29+37+digits
expectedResult: osmomath.MustNewDecFromStr("1.999999999999999999986137056388801094"),

errTolerance: osmoutils.ErrTolerance{
AdditiveTolerance: minDecTolerance,
MultiplicativeTolerance: minDecTolerance,
RoundingDir: osmomath.RoundDown,
},
},
"exp2(0.0000...1)": {
exponent: osmomath.ZeroDec().Add(smallerValue),
// https://www.wolframalpha.com/input?i=2%5E0.00000000000000000001+37+digits
expectedResult: osmomath.MustNewDecFromStr("1.000000000000000000006931471805599453"),

errTolerance: osmoutils.ErrTolerance{
AdditiveTolerance: minDecTolerance,
MultiplicativeTolerance: minDecTolerance,
RoundingDir: osmomath.RoundUnconstrained, // TODO: confirm if this is acceptable.
},
},
"exp2(0.3334567)": {
exponent: osmomath.MustNewDecFromStr("0.3334567"),
// https://www.wolframalpha.com/input?i=2%5E0.3334567+37+digits
expectedResult: osmomath.MustNewDecFromStr("1.260028791934303989065848870753742298"),

errTolerance: osmoutils.ErrTolerance{
AdditiveTolerance: sdk.MustNewDecFromStr("0.00000000000000007"),
MultiplicativeTolerance: minDecTolerance.Mul(sdk.NewDec(10)),
RoundingDir: osmomath.RoundDown,
},
},
"exp2(0.84864288)": {
exponent: osmomath.MustNewDecFromStr("0.84864288"),
// https://www.wolframalpha.com/input?i=2%5E0.84864288+37+digits
expectedResult: osmomath.MustNewDecFromStr("1.800806138872630518880998772777747572"),

errTolerance: osmoutils.ErrTolerance{
AdditiveTolerance: sdk.MustNewDecFromStr("0.00000000000000002"),
MultiplicativeTolerance: minDecTolerance.Mul(sdk.NewDec(10)),
RoundingDir: osmomath.RoundUnconstrained, // TODO: confirm if this is acceptable.
},
},
"exp2(0.999999999999999999999999999999999956)": {
exponent: osmomath.MustNewDecFromStr("0.999999999999999999999999999999999956"),
// https://www.wolframalpha.com/input?i=2%5E0.999999999999999999999999999999999956+37+digits
expectedResult: osmomath.MustNewDecFromStr("1.999999999999999999999999999999999939"),

errTolerance: osmoutils.ErrTolerance{
AdditiveTolerance: minDecTolerance,
MultiplicativeTolerance: minDecTolerance,
RoundingDir: osmomath.RoundDown,
},
},
// out of bounds.
"exponent < 0 - panic": {
exponent: osmomath.ZeroDec().Sub(smallValue),
expectPanic: true,
},
"exponent > 1 - panic": {
exponent: osmomath.OneDec().Add(smallValue),
expectPanic: true,
},
}

for name, tc := range tests {
tc := tc
t.Run(name, func(t *testing.T) {
osmoassert.ConditionalPanic(t, tc.expectPanic, func() {
// System under test.
result := osmomath.Exp2ChebyshevRationalApprox(tc.exponent)

// Reuse the same test cases for exp2 that is a wrapper around Exp2ChebyshevRationalApprox.
// This is done to reduce boilerplate from duplicating test cases.
resultExp2 := osmomath.Exp2(tc.exponent)
require.Equal(t, result, resultExp2)

require.Equal(t, 0, tc.errTolerance.CompareBigDec(tc.expectedResult, result))
})
})
}
}

func TestExp2(t *testing.T) {
tests := map[string]struct {
exponent osmomath.BigDec
expectedResult osmomath.BigDec
errTolerance osmoutils.ErrTolerance
expectPanic bool
}{
"exp2(28.5)": {
exponent: osmomath.MustNewDecFromStr("28.5"),
// https://www.wolframalpha.com/input?i=2%5E%2828.5%29+45+digits
expectedResult: osmomath.MustNewDecFromStr("379625062.497006211556423566253288543343173698"),

errTolerance: osmoutils.ErrTolerance{
AdditiveTolerance: minDecTolerance,
MultiplicativeTolerance: minDecTolerance,
RoundingDir: osmomath.RoundUnconstrained, // TODO: confirm if this is acceptable.
p0mvn marked this conversation as resolved.
Show resolved Hide resolved
},
},
"exp2(63.84864288)": {
exponent: osmomath.MustNewDecFromStr("63.84864288"),
// https://www.wolframalpha.com/input?i=2%5E%2863.84864288%29+56+digits
expectedResult: osmomath.MustNewDecFromStr("16609504985074238416.013387053450559984846024066925604094"),

errTolerance: osmoutils.ErrTolerance{
AdditiveTolerance: sdk.MustNewDecFromStr("0.00042"),
MultiplicativeTolerance: minDecTolerance,
RoundingDir: osmomath.RoundUnconstrained, // TODO: confirm if this is acceptable.
},
},
"exp2(64.5)": {
exponent: osmomath.MustNewDecFromStr("64.5"),
// https://www.wolframalpha.com/input?i=2%5E%2864.5%29+56+digits
expectedResult: osmomath.MustNewDecFromStr("26087635650665564424.699143612505016737766552579185717157"),

errTolerance: osmoutils.ErrTolerance{
AdditiveTolerance: sdk.MustNewDecFromStr("0.000000000000000008"),
MultiplicativeTolerance: minDecTolerance,
RoundingDir: osmomath.RoundUnconstrained, // TODO: confirm if this is acceptable.
},
},
"exp2(80.5)": {
exponent: osmomath.MustNewDecFromStr("80.5"),
// https://www.wolframalpha.com/input?i=2%5E%2880.5%29+61+digits
expectedResult: osmomath.MustNewDecFromStr("1709679290002018430137083.075789128776926268789829515159631571"),

errTolerance: osmoutils.ErrTolerance{
AdditiveTolerance: sdk.MustNewDecFromStr("0.0000000000006"),
MultiplicativeTolerance: minDecTolerance,
RoundingDir: osmomath.RoundUnconstrained, // TODO: confirm if this is acceptable.
},
},
"exp2(100.5)": {
exponent: osmomath.MustNewDecFromStr("100.5"),
// https://www.wolframalpha.com/input?i=2%5E%28100.5%29+67+digits
expectedResult: osmomath.MustNewDecFromStr("1792728671193156477399422023278.661496394239222564273688025833797661"),

errTolerance: osmoutils.ErrTolerance{
AdditiveTolerance: sdk.MustNewDecFromStr("0.0000006"),
MultiplicativeTolerance: minDecTolerance,
RoundingDir: osmomath.RoundUnconstrained, // TODO: confirm if this is acceptable.
},
},
"exp2(128.5)": {
exponent: osmomath.MustNewDecFromStr("128.5"),
// https://www.wolframalpha.com/input?i=2%5E%28128.5%29+75+digits
expectedResult: osmomath.MustNewDecFromStr("481231938336009023090067544955250113854.229961482126296754016435255422777776"),

errTolerance: osmoutils.ErrTolerance{
AdditiveTolerance: sdk.MustNewDecFromStr("146.5"),
MultiplicativeTolerance: minDecTolerance,
RoundingDir: osmomath.RoundUnconstrained, // TODO: confirm if this is acceptable.
},
},
"exp2(127.999999999999999999999999999999999999)": {
exponent: osmomath.MustNewDecFromStr("127.999999999999999999999999999999999999"),
// https://www.wolframalpha.com/input?i=2%5E%28127.999999999999999999999999999999999999%29+75+digits
expectedResult: osmomath.MustNewDecFromStr("340282366920938463463374607431768211220.134236774486705862055857235845515682"),

errTolerance: osmoutils.ErrTolerance{
AdditiveTolerance: sdk.MustNewDecFromStr("15044647266406936"),
MultiplicativeTolerance: minDecTolerance,
RoundingDir: osmomath.RoundDown,
},
},
"exp2(127.84864288)": {
exponent: osmomath.MustNewDecFromStr("127.84864288"),
// https://www.wolframalpha.com/input?i=2%5E%28127.84864288%29+75+digits
expectedResult: osmomath.MustNewDecFromStr("306391287650667462068703337664945630660.398687487527674545778353588077174571"),

errTolerance: osmoutils.ErrTolerance{
AdditiveTolerance: sdk.MustNewDecFromStr("7707157415597963"),
MultiplicativeTolerance: minDecTolerance,
RoundingDir: osmomath.RoundUnconstrained,
},
},
"panic, too large - positive": {
exponent: osmomath.MaxSupportedExponent.Add(osmomath.OneDec()),
expectPanic: true,
},
"panic - negative exponent": {
exponent: osmomath.OneDec().Neg(),
expectPanic: true,
},
"at exponent boundary - positive": {
exponent: osmomath.MaxSupportedExponent,
// https://www.wolframalpha.com/input?i=2%5E%282%5E9%29
expectedResult: osmomath.MustNewDecFromStr("13407807929942597099574024998205846127479365820592393377723561443721764030073546976801874298166903427690031858186486050853753882811946569946433649006084096"),

errTolerance: osmoutils.ErrTolerance{
AdditiveTolerance: minDecTolerance,
MultiplicativeTolerance: minDecTolerance,
RoundingDir: osmomath.RoundDown,
},
},
}

for name, tc := range tests {
tc := tc
t.Run(name, func(t *testing.T) {
osmoassert.ConditionalPanic(t, tc.expectPanic, func() {

// System under test.
result := osmomath.Exp2(tc.exponent)

require.Equal(t, 0, tc.errTolerance.CompareBigDec(tc.expectedResult, result))
})
})
}
}
Loading