Repository: osmosis-labs/osmosis
Pull request #3708 (2022-12-14 00:43:45)

feat(osmomath): Exp2 function
<!-- < < < < < < < < < < < < < < < < < < < < < < < < < < < < < < < < < ☺
v                               ✰  Thanks for creating a PR! ✰    
v    Before smashing the submit button please review the checkboxes.
v    If a checkbox is n/a - please still include it but + a little note why
v    If your PR doesn't close an issue, that's OK!  Just remove the Closes: #XXX line!
☺ > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > >  -->

Closes: #XXX

## What is the purpose of the change

This PR adds an `Exp2` function that computes $$2^{decimal}$$

The following property is used:
$$2^{decimal} = 2^{integer} * 2^{fractional}$$

The integer exponent result is computed using the existing `PowerInteger` function.

The fractional exponent result is computed using the 13-parameter Chebyshev rational approximation. The numerator coefficients are rounded down and the denominator coefficients are rounded up at the precision end. However, this does not guarantee that the approximated result is always smaller than the actual results (as determined by WolframAlpha).

The Python approximation script is modified to utilize the [production coefficients ](https://github.com/osmosis-labs/osmosis/blob/0cc8250f2795860ef2fe41896a0dc44ff68005dd/scripts/approximations/main.py#L194-L212) for estimating relative errors. Please see the results:
![image](https://user-images.githubusercontent.com/34196718/207486684-f3656df7-0051-4d34-b6bf-95bc483abc6f.png)

Additionally, `ErrTolerance` was refactored to use `sdk.Dec` for `AdditiveTolerance` for more accuracy in tests.

In a future PR, `Power` function with a decimal exponent can be added. This function will use `exp2` by utilizing the following property:
`x^y = exp2(y * log_2{x})`

## Testing and Verifying

Tests with links to WolframAlpha were added. The absolute and the relative errors are bounded.

## Open Questions

- Is it acceptable that we cannot guarantee rounding behavior?

## Documentation and Release Note

  - Does this pull request introduce a new feature or user-facing behavior changes? no
  - Is a relevant changelog entry added to the `Unreleased` section in `CHANGELOG.md`? yes
  - How is the feature or change documented? godocs and approximation scripts

nicolaslara: > Is it acceptable that we cannot guarantee rounding behavior?

I think as long as the error is within the an acceptable range (from the chart $[10^{-22}, 10^{-27}]$), then I think it doesn't matter if it's higher or lower than the actual results
p0mvn: > > Is it acceptable that we cannot guarantee rounding behavior?
> 
> I think as long as the error is within the an acceptable range (from the chart [10−22,10−27]), then I think it doesn't matter if it's higher or lower than the actual results

There are some cases where the error is outside of the chart - `10^17`. For example:
https://github.com/osmosis-labs/osmosis/blob/79beaae72f54159a82c54a7872b6053fc8ffdc5a/osmomath/exp2_test.go#L112

I have not been able to explain why that is the case.

One guess is that approximated results in Go are being compared against WolframAlpha. The chart is constructed from a Python library where both the approximations and the actual results are using the same framework. I'm wondering if different computational platforms introduce additional errors.

Please let me know your thoughts.
p0mvn: This is how the relative error plot looks like (very similar if not identical to the absolute):
![Screenshot from 2022-12-15 14-53-00](https://user-images.githubusercontent.com/34196718/207954194-fb0e497a-38ea-4d7e-b122-e195fb672ca4.png)

From call, the concern was that relative and absolute plots should not look the same.

On a quick review, I do not see a bug in the way it is calculated. However, I'm going to investigate this better and get back with updates.
ValarDragon: Python and actual go code both LGTM, nice job! Still reviewing go tests
ValarDragon: The plot difference LGTM

Were actually seeing what we'd expect to see between the absolute error vs relative error

The additive error on RHS side was larger than LHS. The relative error 'dampens' the RHS relative to the LHS, as exp only increasing on the range
p0mvn: > The plot difference LGTM
> 
> Were actually seeing what we'd expect to see between the absolute error vs relative error
> 
> The additive error on RHS side was larger than LHS. The relative error 'dampens' the RHS relative to the LHS, as exp only increasing on the range

That makes sense. Thanks for confirming
p0mvn: All comments are now addressed. This should be ready for merge once CI passes.

Patch:
diff --git a/CHANGELOG.md b/CHANGELOG.md
index a1e8ebb5a..e713afb38 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -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
 
 ### API breaks
diff --git a/osmomath/exp2.go b/osmomath/exp2.go
new file mode 100644
index 000000000..b7973b3bc
--- /dev/null
+++ b/osmomath/exp2.go
@@ -0,0 +1,101 @@
+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"),
+	}
+
+	// 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.
+// The answer is correct up to a factor of 10^-18.
+// Meaning, result = result * k for k in [1 - 10^(-18), 1 + 10^(-18)]
+// Note: our Python script plots show accuracy up to a factor of 10^22.
+// However, in Go tests we only test up to 10^18. Therefore, this is the guarantee.
+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)
+
+	// Left bit shift is equivalent to multiplying by 2^integerExponent.
+	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.
+// The answer is correct up to a factor of 10^-18.
+// Meaning, result = result * k for k in [1 - 10^(-18), 1 + 10^(-18)]
+// Note: our Python script plots show accuracy up to a factor of 10^22.
+// However, in Go tests we only test up to 10^18. Therefore, this is the guarantee.
+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)
+}
diff --git a/scripts/approximations/approximations.py b/scripts/approximations/approximations.py
index 198496980..1a4345096 100644
--- a/scripts/approximations/approximations.py
+++ b/scripts/approximations/approximations.py
@@ -84,9 +84,9 @@ def compute_max_error(y_approximation, y_actual) -> sp.Float:
             max = sp.Max(max, cur_abs)
     return max
 
-def compute_error_range(y_approximation: list, y_actual: list) -> list[sp.Float]:
+def compute_absolute_error_range(y_approximation: list, y_actual: list) -> list[sp.Float]:
     """ Given an approximated list of y values and actual y values, computes and returns
-    error deltas between them.
+    absolute error between them, computed as | y_approximation[i] - y_actual[i] |.
 
     CONTRACT:
     - y_approximation and y_actual must be the same length
@@ -101,6 +101,29 @@ def compute_error_range(y_approximation: list, y_actual: list) -> list[sp.Float]
         result.append(cur_abs)
     return result
 
+def compute_relative_error_range(y_approximation: list, y_actual: list) -> list[sp.Float]:
+    """ Given an approximated list of y values and actual y values, computes and returns
+    relative error between them, computed as | y_approximation[i] - y_actual[i] | / y_actual[i].
+
+    For y_actual[i] = 0, relative error is defined as 0.
+
+    CONTRACT:
+    - y_approximation and y_actual must be the same length
+    - for every i in range(len(y_approximation)), y_approximation[i] and y_actual[i] must correspond to the
+    same x coordinate
+    """
+    result = []
+    for i in range(len(y_approximation)):
+        if y_actual[i] == 0:
+            result.append(0)
+            continue
+
+        cur_relative_error = sp.Abs(y_approximation[i] - y_actual[i]) / y_actual[i]
+        if cur_relative_error is sp.nan:
+            raise ValueError(F"cur_abs is nan. y_approximation[i] ({y_approximation[i]}) and y_actual[i] ({y_actual[i]})")
+        result.append(cur_relative_error)
+    return result
+
 def equispaced_poly_approx(fn, x_start: sp.Float, x_end: sp.Float, num_terms: int):
     """ Returns the coefficients for an equispaced polynomial between x_start and x_end with num_terms terms.
 
diff --git a/scripts/approximations/main.py b/scripts/approximations/main.py
index b91433855..47875ed7b 100644
--- a/scripts/approximations/main.py
+++ b/scripts/approximations/main.py
@@ -20,10 +20,10 @@ num_parameters = 13
 num_parameters_errors = 30
 
 # number of (x,y) coordinates used to plot the resulting approximation.
-num_points_plot = 100000
+num_points_plot = 10000
 
 # function to approximate
-approximated_fn = lambda x: sp.Pow(sp.E, x)
+approximated_fn = lambda x: sp.Pow(2, x)
 
 # fixed point precision used in Osmosis `osmomath` package.
 osmomath_precision = 36
@@ -48,16 +48,23 @@ shouldPlotErrorRange  = True
 # Plots if true.
 shouldPlotMaxError = True
 
-def plot_error_range(x_coordinates, y_approximation, y_actual):
+def plot_error_range(x_coordinates, y_approximation, y_actual, is_absolute: bool):
     """ Given x coordinates that correspond to approximated y coordinates and actual y coordinates,
-    compute the deltas between y approximated and y actual and plot them in log scale on y.
+    computes the error between y approximated and y actual and plot them in log scale on y.
+
+    If is_absolute, plots the absolute error. Otherwise, plots the relative error.
     """
-    error_deltas = approximations.compute_error_range(y_approximation, y_actual)
+    if is_absolute:
+        error_kind_str = "Absolute"
+        error_deltas = approximations.compute_absolute_error_range(y_approximation, y_actual)
+    else:
+        error_kind_str = "Relative"
+        error_deltas = approximations.compute_relative_error_range(y_approximation, y_actual)
 
     plt.semilogy(x_coordinates, error_deltas)
 
     plt.grid(True)
-    plt.title(f"Chebyshev Rational e^x Errors on [{x_start}, {x_end}]. {num_parameters} params, {num_points_plot} points")
+    plt.title(f"Chebyshev Rational e^x {error_kind_str} Errors on [{x_start}, {x_end}]. {num_parameters} params, {num_points_plot} points")
     plt.show()
 
 # This script does the following:
@@ -171,14 +178,41 @@ def main():
 def exponent_approximation_choice():
     # Equispaced x coordinates to be used for plotting every approximation.
     x_coordinates = approximations.linspace(x_start, x_end, num_points_plot)
-    x_coordinates = [sp.Float(sp.N(coef, osmomath_precision + 1), osmomath_precision + 1) for coef in x_coordinates]
+    x_coordinates = [sp.N(sp.Float(coef, osmomath_precision), n=osmomath_precision) for coef in x_coordinates]
 
-    # Chebyshev Rational Approximation to get the coefficients.
-    coef_numerator, coef_denominator = approximations.chebyshev_rational_approx(approximated_fn, x_start, x_end, num_parameters)
+    print(x_coordinates)
 
-    # Truncate the coefficients to osmomath precision.
-    coef_numerator = [sp.Float(sp.N(coef, osmomath_precision + 1), osmomath_precision + 1) for coef in coef_numerator]
-    coef_denominator = [sp.Float(sp.N(coef, osmomath_precision + 1), osmomath_precision + 1) for coef in coef_denominator]
+    # Chebyshev Rational Approximation to get the coefficients.
+    # coef_numerator, coef_denominator = approximations.chebyshev_rational_approx(approximated_fn, x_start, x_end, num_parameters)
+    # coef_numerator = [sp.N(coef, osmomath_precision + 2) for coef in coef_numerator]
+    # coef_denominator = [sp.N(coef, osmomath_precision + 2) for coef in coef_denominator]
+
+    # Hard code and round up numerator coefficientst that are to be used in production
+    # Hard code and round down numerator coefficientst that are to be used in production
+    # Both of these are calculated us=ing the above commented out code.
+
+    coef_numerator = [
+        sp.N(sp.Float("1.000000000000000000000044212244679434", osmomath_precision), n=osmomath_precision),
+        sp.N(sp.Float("0.352032455817400196452603772766844426", osmomath_precision), n=osmomath_precision),
+        sp.N(sp.Float("0.056507868883666405413116800969512484", osmomath_precision), n=osmomath_precision),
+        sp.N(sp.Float("0.005343900728213034434757419480319916", osmomath_precision), n=osmomath_precision),
+        sp.N(sp.Float("0.000317708814342353603087543715930732", osmomath_precision), n=osmomath_precision),
+        sp.N(sp.Float("0.000011429747507407623028722262874632", osmomath_precision), n=osmomath_precision),
+        sp.N(sp.Float("0.000000198381965651614980168744540366", osmomath_precision), n=osmomath_precision),
+    ]
+
+    coef_denominator = [
+        sp.N(sp.Float("1.0000000000000000000000000000000000000", osmomath_precision), n=osmomath_precision),
+        sp.N(sp.Float("-0.341114724742545112949699755780593311", osmomath_precision), n=osmomath_precision),
+        sp.N(sp.Float("0.052724071627342653404436933178482287", osmomath_precision), n=osmomath_precision),
+        sp.N(sp.Float("-0.004760950735524957576233524801866342", osmomath_precision), n=osmomath_precision),
+        sp.N(sp.Float("0.000267168475410566529819971616894193", osmomath_precision), n=osmomath_precision),
+        sp.N(sp.Float("-0.000008923715368802211181557353097439", osmomath_precision), n=osmomath_precision),
+        sp.N(sp.Float("0.000000140277233177373698516010555916", osmomath_precision), n=osmomath_precision),
+    ]
+
+    print(coef_numerator)
+    print(coef_denominator)
 
     # Evaluate approximation.
     y_chebyshev_rational = rational.evaluate(x_coordinates, coef_numerator, coef_denominator)
@@ -186,7 +220,7 @@ def exponent_approximation_choice():
     # Compute Actual Values
     y_actual = approximations.get_y_actual(approximated_fn, x_coordinates)
 
-    plot_error_range(x_coordinates, y_chebyshev_rational, y_actual)
+    plot_error_range(x_coordinates, y_chebyshev_rational, y_actual, True)
 
 if __name__ == "__main__":
     # Uncomment to run the main script.


Test patch:
diff --git a/osmomath/exp2_test.go b/osmomath/exp2_test.go
new file mode 100644
index 000000000..0955251dd
--- /dev/null
+++ b/osmomath/exp2_test.go
@@ -0,0 +1,301 @@
+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"
+)
+
+var (
+	// minDecTolerance minimum tolerance for sdk.Dec, given its precision of 18.
+	minDecTolerance = sdk.MustNewDecFromStr("0.000000000000000001")
+)
+
+func TestExp2ChebyshevRationalApprox(t *testing.T) {
+	// These values are used to test the approximated results close
+	// to 0 and 1 boundaries.
+	// With other types of approximations, there is a high likelyhood
+	// of larger errors clsoer to the boundaries. This is known as Runge's phenomenon.
+	// https://en.wikipedia.org/wiki/Runge%27s_phenomenon
+	//
+	// Chebyshev approximation should be able to handle this better.
+	// Tests at the boundaries help to validate there is no Runge's phenomenon.
+	smallValue := osmomath.MustNewDecFromStr("0.00001")
+	smallerValue := osmomath.MustNewDecFromStr("0.00000000000000000001")
+
+	tests := map[string]struct {
+		exponent       osmomath.BigDec
+		expectedResult osmomath.BigDec
+		errTolerance   osmomath.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: osmomath.ErrTolerance{
+				AdditiveTolerance:       minDecTolerance,
+				MultiplicativeTolerance: minDecTolerance,
+				RoundingDir:             osmomath.RoundDown,
+			},
+		},
+		"exp2(0)": {
+			exponent:       osmomath.ZeroDec(),
+			expectedResult: osmomath.OneDec(),
+
+			errTolerance: osmomath.ErrTolerance{
+				AdditiveTolerance:       sdk.ZeroDec(),
+				MultiplicativeTolerance: sdk.ZeroDec(),
+				RoundingDir:             osmomath.RoundDown,
+			},
+		},
+		"exp2(1)": {
+			exponent:       osmomath.OneDec(),
+			expectedResult: osmomath.MustNewDecFromStr("2"),
+
+			errTolerance: osmomath.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: osmomath.ErrTolerance{
+				AdditiveTolerance:       minDecTolerance,
+				MultiplicativeTolerance: minDecTolerance,
+				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: osmomath.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: osmomath.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: osmomath.ErrTolerance{
+				AdditiveTolerance:       minDecTolerance,
+				MultiplicativeTolerance: minDecTolerance,
+				RoundingDir:             osmomath.RoundUnconstrained,
+			},
+		},
+		"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: osmomath.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: osmomath.ErrTolerance{
+				AdditiveTolerance:       sdk.MustNewDecFromStr("0.00000000000000002"),
+				MultiplicativeTolerance: minDecTolerance.Mul(sdk.NewDec(10)),
+				RoundingDir:             osmomath.RoundUnconstrained,
+			},
+		},
+		"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: osmomath.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   osmomath.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: osmomath.ErrTolerance{
+				AdditiveTolerance:       minDecTolerance,
+				MultiplicativeTolerance: minDecTolerance,
+				RoundingDir:             osmomath.RoundUnconstrained,
+			},
+		},
+		"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: osmomath.ErrTolerance{
+				AdditiveTolerance:       sdk.MustNewDecFromStr("0.00042"),
+				MultiplicativeTolerance: minDecTolerance,
+				RoundingDir:             osmomath.RoundUnconstrained,
+			},
+		},
+		"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: osmomath.ErrTolerance{
+				AdditiveTolerance:       sdk.MustNewDecFromStr("0.000000000000000008"),
+				MultiplicativeTolerance: minDecTolerance,
+				RoundingDir:             osmomath.RoundUnconstrained,
+			},
+		},
+		"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: osmomath.ErrTolerance{
+				AdditiveTolerance:       sdk.MustNewDecFromStr("0.0000000000006"),
+				MultiplicativeTolerance: minDecTolerance,
+				RoundingDir:             osmomath.RoundUnconstrained,
+			},
+		},
+		"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: osmomath.ErrTolerance{
+				AdditiveTolerance:       sdk.MustNewDecFromStr("0.0000006"),
+				MultiplicativeTolerance: minDecTolerance,
+				RoundingDir:             osmomath.RoundUnconstrained,
+			},
+		},
+		"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: osmomath.ErrTolerance{
+				AdditiveTolerance:       sdk.MustNewDecFromStr("146.5"),
+				MultiplicativeTolerance: minDecTolerance,
+				RoundingDir:             osmomath.RoundUnconstrained,
+			},
+		},
+		"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: osmomath.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: osmomath.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: osmomath.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))
+			})
+		})
+	}
+}
diff --git a/osmomath/export_test.go b/osmomath/export_test.go
index bdb6c71f9..d723e1060 100644
--- a/osmomath/export_test.go
+++ b/osmomath/export_test.go
@@ -8,8 +8,9 @@ import (
 )
 
 var (
-	EulersNumber = eulersNumber
-	TwoBigDec    = twoBigDec
+	MaxSupportedExponent = maxSupportedExponent
+	EulersNumber         = eulersNumber
+	TwoBigDec            = twoBigDec
 )
 
 // 2^128 - 1, needs to be the same as gammtypes.MaxSpotPrice
@@ -28,3 +29,7 @@ func ConditionalPanic(t *testing.T, expectPanic bool, sut func()) {
 	}
 	require.NotPanics(t, sut)
 }
+
+func Exp2ChebyshevRationalApprox(exponent BigDec) BigDec {
+	return exp2ChebyshevRationalApprox(exponent)
+}