From a1ef64da45f2e0fef49e52ddf530bb38ffe8bdbc Mon Sep 17 00:00:00 2001 From: DrunkRandomWalker Date: Thu, 21 Sep 2023 10:15:51 -0400 Subject: [PATCH] minor improvements --- packages/injective-math/src/fp_decimal/exp.rs | 198 ++++++++++++++---- packages/injective-math/src/fp_decimal/log.rs | 55 ++--- 2 files changed, 190 insertions(+), 63 deletions(-) diff --git a/packages/injective-math/src/fp_decimal/exp.rs b/packages/injective-math/src/fp_decimal/exp.rs index 55800d2e..31e7736b 100644 --- a/packages/injective-math/src/fp_decimal/exp.rs +++ b/packages/injective-math/src/fp_decimal/exp.rs @@ -360,24 +360,132 @@ impl FPDecimal { } } - fn inner(mut a: FPDecimal, mut b: FPDecimal) -> Result { + fn inner(mut a: FPDecimal, mut exponent: FPDecimal) -> Result { // a^b // 14 terms taylor expansion provides a good enough approximation let n_terms = 13u128; - match b.cmp(&FPDecimal::ZERO) { + match exponent.cmp(&FPDecimal::ZERO) { Ordering::Equal => Ok(FPDecimal::one()), Ordering::Less => { - a = FPDecimal::ONE / a; - b = -b; - match b.cmp(&(FPDecimal::ONE)) { - Ordering::Equal => Ok(a), - Ordering::Less => Ok(FPDecimal::_exp_taylor_expansion(a, b, n_terms)), + exponent = -exponent; + match exponent.cmp(&(FPDecimal::ONE)) { + Ordering::Equal => Ok(FPDecimal::ONE / a), + Ordering::Less => { + // NOTE: only accurate for 1,3,5,7,11, and combinations of these numbers + if a.log2().is_some() { + if ((FPDecimal::reciprocal(exponent) % FPDecimal::TWO).int() - FPDecimal::ONE).abs() + <= FPDecimal::must_from_str("0.000001") + { + let mut tmp_b = FPDecimal::reciprocal(exponent).int(); + while tmp_b > FPDecimal::ONE { + a /= FPDecimal::TWO; + tmp_b -= FPDecimal::ONE; + } + return Ok(FPDecimal::ONE / a); + }; + if FPDecimal::reciprocal(exponent) % FPDecimal::TWO == FPDecimal::ZERO { + let mut tmp_b = FPDecimal::reciprocal(exponent).int(); + while tmp_b > FPDecimal::ONE { + a = a.sqrt(); + tmp_b /= FPDecimal::TWO; + } + return Ok(FPDecimal::ONE / a); + }; + } + + if a.log3().is_some() { + if ((FPDecimal::reciprocal(exponent) % FPDecimal::TWO).int() - FPDecimal::ONE).abs() + <= FPDecimal::must_from_str("0.000001") + { + let mut tmp_b = FPDecimal::reciprocal(exponent).int(); + while tmp_b > FPDecimal::ONE { + a /= FPDecimal::THREE; + tmp_b -= FPDecimal::ONE; + } + return Ok(FPDecimal::ONE / a); + }; + if FPDecimal::reciprocal(exponent) % FPDecimal::TWO == FPDecimal::ZERO { + let mut tmp_b = FPDecimal::reciprocal(exponent).int(); + while tmp_b > FPDecimal::ONE { + a = a.sqrt(); + tmp_b /= FPDecimal::TWO; + } + return Ok(FPDecimal::ONE / a); + }; + } + + if a.log5().is_some() { + println!("here"); + if ((FPDecimal::reciprocal(exponent) % FPDecimal::TWO).int() - FPDecimal::ONE).abs() + <= FPDecimal::must_from_str("0.000001") + { + let mut tmp_b = FPDecimal::reciprocal(exponent).int(); + while tmp_b > FPDecimal::ONE { + a /= FPDecimal::FIVE; + tmp_b -= FPDecimal::ONE; + } + return Ok(FPDecimal::ONE / a); + }; + if FPDecimal::reciprocal(exponent) % FPDecimal::TWO == FPDecimal::ZERO { + let mut tmp_b = FPDecimal::reciprocal(exponent).int(); + while tmp_b > FPDecimal::ONE { + a = a.sqrt(); + tmp_b /= FPDecimal::TWO; + } + return Ok(FPDecimal::ONE / a); + }; + } + + if a.log7().is_some() { + if ((FPDecimal::reciprocal(exponent) % FPDecimal::TWO).int() - FPDecimal::ONE).abs() + <= FPDecimal::must_from_str("0.000001") + { + let mut tmp_b = FPDecimal::reciprocal(exponent).int(); + while tmp_b > FPDecimal::ONE { + a /= FPDecimal::SEVEN; + tmp_b -= FPDecimal::ONE; + } + return Ok(FPDecimal::ONE / a); + }; + if FPDecimal::reciprocal(exponent) % FPDecimal::TWO == FPDecimal::ZERO { + let mut tmp_b = FPDecimal::reciprocal(exponent).int(); + while tmp_b > FPDecimal::ONE { + a = a.sqrt(); + tmp_b /= FPDecimal::TWO; + } + return Ok(FPDecimal::ONE / a); + }; + } + + if a.log11().is_some() { + if ((FPDecimal::reciprocal(exponent) % FPDecimal::TWO).int() - FPDecimal::ONE).abs() + <= FPDecimal::must_from_str("0.000001") + { + let mut tmp_b = FPDecimal::reciprocal(exponent).int(); + while tmp_b > FPDecimal::ONE { + a /= FPDecimal::from(11u128); + tmp_b -= FPDecimal::ONE; + } + return Ok(FPDecimal::ONE / a); + }; + if FPDecimal::reciprocal(exponent) % FPDecimal::TWO == FPDecimal::ZERO { + let mut tmp_b = FPDecimal::reciprocal(exponent).int(); + while tmp_b > FPDecimal::ONE { + a = a.sqrt(); + tmp_b /= FPDecimal::TWO; + } + return Ok(FPDecimal::ONE / a); + }; + } + + Ok(FPDecimal::_exp_taylor_expansion(FPDecimal::ONE / a, exponent, n_terms)) + } Ordering::Greater => { - let mut int_b = b.int(); - let rem_b = b - int_b; + let mut int_b = exponent.int(); + let rem_b = exponent - int_b; let mut float_exp = FPDecimal::ONE; if rem_b != FPDecimal::ZERO { - float_exp = FPDecimal::_exp_taylor_expansion(a, rem_b, n_terms); + float_exp = FPDecimal::_exp_taylor_expansion(FPDecimal::ONE / a, rem_b, n_terms); } let mut tmp_a = FPDecimal::ONE; while int_b > FPDecimal::one() { @@ -389,28 +497,29 @@ impl FPDecimal { int_b /= FPDecimal::TWO; } a *= tmp_a; - a *= float_exp; - Ok(a) + // a *= float_exp; + Ok(FPDecimal::ONE / a * float_exp) } } } - Ordering::Greater => match b.cmp(&FPDecimal::ONE) { + Ordering::Greater => match exponent.cmp(&FPDecimal::ONE) { Ordering::Equal => Ok(a), Ordering::Less => { // taylor expansion approximation of exponentation compuation with float number exponent // NOTE: only accurate for 1,3,5,7,11, and combinations of these numbers if a.log2().is_some() { - if ((FPDecimal::reciprocal(b) % FPDecimal::TWO).int() - FPDecimal::ONE).abs() <= FPDecimal::must_from_str("0.000001") + if ((FPDecimal::reciprocal(exponent) % FPDecimal::TWO).int() - FPDecimal::ONE).abs() + <= FPDecimal::must_from_str("0.000001") { - let mut tmp_b = FPDecimal::reciprocal(b).int(); + let mut tmp_b = FPDecimal::reciprocal(exponent).int(); while tmp_b > FPDecimal::ONE { a /= FPDecimal::TWO; tmp_b -= FPDecimal::ONE; } return Ok(a); }; - if FPDecimal::reciprocal(b) % FPDecimal::TWO == FPDecimal::ZERO { - let mut tmp_b = FPDecimal::reciprocal(b).int(); + if FPDecimal::reciprocal(exponent) % FPDecimal::TWO == FPDecimal::ZERO { + let mut tmp_b = FPDecimal::reciprocal(exponent).int(); while tmp_b > FPDecimal::ONE { a = a.sqrt(); tmp_b /= FPDecimal::TWO; @@ -420,17 +529,18 @@ impl FPDecimal { } if a.log3().is_some() { - if ((FPDecimal::reciprocal(b) % FPDecimal::TWO).int() - FPDecimal::ONE).abs() <= FPDecimal::must_from_str("0.000001") + if ((FPDecimal::reciprocal(exponent) % FPDecimal::TWO).int() - FPDecimal::ONE).abs() + <= FPDecimal::must_from_str("0.000001") { - let mut tmp_b = FPDecimal::reciprocal(b).int(); + let mut tmp_b = FPDecimal::reciprocal(exponent).int(); while tmp_b > FPDecimal::ONE { a /= FPDecimal::THREE; tmp_b -= FPDecimal::ONE; } return Ok(a); }; - if FPDecimal::reciprocal(b) % FPDecimal::TWO == FPDecimal::ZERO { - let mut tmp_b = FPDecimal::reciprocal(b).int(); + if FPDecimal::reciprocal(exponent) % FPDecimal::TWO == FPDecimal::ZERO { + let mut tmp_b = FPDecimal::reciprocal(exponent).int(); while tmp_b > FPDecimal::ONE { a = a.sqrt(); tmp_b /= FPDecimal::TWO; @@ -441,17 +551,18 @@ impl FPDecimal { if a.log5().is_some() { println!("here"); - if ((FPDecimal::reciprocal(b) % FPDecimal::TWO).int() - FPDecimal::ONE).abs() <= FPDecimal::must_from_str("0.000001") + if ((FPDecimal::reciprocal(exponent) % FPDecimal::TWO).int() - FPDecimal::ONE).abs() + <= FPDecimal::must_from_str("0.000001") { - let mut tmp_b = FPDecimal::reciprocal(b).int(); + let mut tmp_b = FPDecimal::reciprocal(exponent).int(); while tmp_b > FPDecimal::ONE { a /= FPDecimal::FIVE; tmp_b -= FPDecimal::ONE; } return Ok(a); }; - if FPDecimal::reciprocal(b) % FPDecimal::TWO == FPDecimal::ZERO { - let mut tmp_b = FPDecimal::reciprocal(b).int(); + if FPDecimal::reciprocal(exponent) % FPDecimal::TWO == FPDecimal::ZERO { + let mut tmp_b = FPDecimal::reciprocal(exponent).int(); while tmp_b > FPDecimal::ONE { a = a.sqrt(); tmp_b /= FPDecimal::TWO; @@ -461,17 +572,18 @@ impl FPDecimal { } if a.log7().is_some() { - if ((FPDecimal::reciprocal(b) % FPDecimal::TWO).int() - FPDecimal::ONE).abs() <= FPDecimal::must_from_str("0.000001") + if ((FPDecimal::reciprocal(exponent) % FPDecimal::TWO).int() - FPDecimal::ONE).abs() + <= FPDecimal::must_from_str("0.000001") { - let mut tmp_b = FPDecimal::reciprocal(b).int(); + let mut tmp_b = FPDecimal::reciprocal(exponent).int(); while tmp_b > FPDecimal::ONE { a /= FPDecimal::SEVEN; tmp_b -= FPDecimal::ONE; } return Ok(a); }; - if FPDecimal::reciprocal(b) % FPDecimal::TWO == FPDecimal::ZERO { - let mut tmp_b = FPDecimal::reciprocal(b).int(); + if FPDecimal::reciprocal(exponent) % FPDecimal::TWO == FPDecimal::ZERO { + let mut tmp_b = FPDecimal::reciprocal(exponent).int(); while tmp_b > FPDecimal::ONE { a = a.sqrt(); tmp_b /= FPDecimal::TWO; @@ -481,17 +593,18 @@ impl FPDecimal { } if a.log11().is_some() { - if ((FPDecimal::reciprocal(b) % FPDecimal::TWO).int() - FPDecimal::ONE).abs() <= FPDecimal::must_from_str("0.000001") + if ((FPDecimal::reciprocal(exponent) % FPDecimal::TWO).int() - FPDecimal::ONE).abs() + <= FPDecimal::must_from_str("0.000001") { - let mut tmp_b = FPDecimal::reciprocal(b).int(); + let mut tmp_b = FPDecimal::reciprocal(exponent).int(); while tmp_b > FPDecimal::ONE { a /= FPDecimal::from(11u128); tmp_b -= FPDecimal::ONE; } return Ok(a); }; - if FPDecimal::reciprocal(b) % FPDecimal::TWO == FPDecimal::ZERO { - let mut tmp_b = FPDecimal::reciprocal(b).int(); + if FPDecimal::reciprocal(exponent) % FPDecimal::TWO == FPDecimal::ZERO { + let mut tmp_b = FPDecimal::reciprocal(exponent).int(); while tmp_b > FPDecimal::ONE { a = a.sqrt(); tmp_b /= FPDecimal::TWO; @@ -500,12 +613,12 @@ impl FPDecimal { }; } - Ok(FPDecimal::_exp_taylor_expansion(a, b, n_terms)) + Ok(FPDecimal::_exp_taylor_expansion(a, exponent, n_terms)) } Ordering::Greater => { - let mut int_b = b.int(); - let rem_b = b - int_b; + let mut int_b = exponent.int(); + let rem_b = exponent - int_b; let mut float_exp = FPDecimal::ONE; if rem_b != FPDecimal::ZERO { float_exp = FPDecimal::_exp_taylor_expansion(a, rem_b, n_terms); @@ -819,7 +932,7 @@ mod tests { let base = FPDecimal::E; assert_eq!( base.pow(FPDecimal::from_str("-3").unwrap()), - FPDecimal::from_str("0.049787068367863942").unwrap() + FPDecimal::from_str("0.049787068367863943").unwrap() ); } @@ -888,7 +1001,7 @@ mod tests { let exponent = FPDecimal::must_from_str("-3.7"); let result = FPDecimal::checked_pow(base, exponent).unwrap(); - assert_eq!(result, FPDecimal::must_from_str("0.045878267230507923")); + assert_eq!(result, FPDecimal::must_from_str("0.045878267230507924")); } #[test] @@ -908,4 +1021,13 @@ mod tests { let result = FPDecimal::checked_pow(base, exponent).unwrap(); assert_eq!(result, FPDecimal::must_from_str("0.716652904902854162")); } + + #[test] + fn test_1_over_16_pow_neg_0_5() { + let base = FPDecimal::ONE / FPDecimal::from(16u128); + let exponent = FPDecimal::must_from_str("-0.5"); + + let result = FPDecimal::checked_pow(base, exponent).unwrap(); + assert_eq!(result, FPDecimal::FOUR); + } } diff --git a/packages/injective-math/src/fp_decimal/log.rs b/packages/injective-math/src/fp_decimal/log.rs index 7c1678c1..bba31cc2 100644 --- a/packages/injective-math/src/fp_decimal/log.rs +++ b/packages/injective-math/src/fp_decimal/log.rs @@ -8,19 +8,15 @@ impl FPDecimal { } if let Some(value) = self.log3() { return Some((value, 3u128)); - // return Some(value); } if let Some(value) = self.log5() { return Some((value, 5u128)); - // return Some(value); } if let Some(value) = self.log7() { return Some((value, 7u128)); - // return Some(value); } if let Some(value) = self.log11() { return Some((value, 11u128)); - // return Some(value); } None } @@ -59,9 +55,36 @@ impl FPDecimal { if self == FPDecimal::from(1024u128) { return Some(FPDecimal::TEN); } - // Some(self.log(FPDecimal::TWO)) - // Some(self.ln() / FPDecimal::TWO.ln()) - // Some(self.log(FPDecimal::TWO)) + if self == FPDecimal::ONE / FPDecimal::TWO { + return Some(-FPDecimal::ONE); + } + if self == FPDecimal::ONE / FPDecimal::FOUR { + return Some(-FPDecimal::TWO); + } + if self == FPDecimal::ONE / FPDecimal::EIGHT { + return Some(-FPDecimal::THREE); + } + if self == FPDecimal::ONE / FPDecimal::from(16u128) { + return Some(-FPDecimal::FOUR); + } + if self == FPDecimal::ONE / FPDecimal::from(32u128) { + return Some(-FPDecimal::FIVE); + } + if self == FPDecimal::ONE / FPDecimal::from(64u128) { + return Some(-FPDecimal::SIX); + } + if self == FPDecimal::ONE / FPDecimal::from(128u128) { + return Some(-FPDecimal::SEVEN); + } + if self == FPDecimal::ONE / FPDecimal::from(256u128) { + return Some(-FPDecimal::EIGHT); + } + if self == FPDecimal::ONE / FPDecimal::from(512u128) { + return Some(-FPDecimal::NINE); + } + if self == FPDecimal::ONE / FPDecimal::from(1024u128) { + return Some(-FPDecimal::TEN); + } None } @@ -130,9 +153,6 @@ impl FPDecimal { if self == FPDecimal::ONE / FPDecimal::from(59049u128) { return Some(-FPDecimal::TEN); } - - // Some(self.ln() / FPDecimal::THREE.ln()) - // Some(self.log(FPDecimal::THREE)) None } @@ -200,9 +220,6 @@ impl FPDecimal { if self == FPDecimal::ONE / FPDecimal::from(9765625u128) { return Some(-FPDecimal::TEN); } - - // Some(self.ln() / FPDecimal::FIVE.ln()) - // Some(self.log(FPDecimal::FIVE)) None } @@ -271,8 +288,6 @@ impl FPDecimal { if self == FPDecimal::ONE / FPDecimal::from(282475249u128) { return Some(-FPDecimal::TEN); } - // Some(self.log(FPDecimal::SEVEN)) - // Some(self.ln() / FPDecimal::SEVEN.ln()) None } @@ -341,9 +356,6 @@ impl FPDecimal { if self == FPDecimal::ONE / FPDecimal::from(25937424601u128) { return Some(-FPDecimal::TEN); } - - // Some(self.log(FPDecimal::from(11u128))) - // Some(self.ln() / FPDecimal::from(11u128).ln()) None } @@ -360,7 +372,6 @@ impl FPDecimal { a.ln() / base.ln() } - // } /// natural logarithm #[allow(clippy::many_single_char_names)] @@ -387,14 +398,9 @@ impl FPDecimal { } if v == FPDecimal::ONE.num { return r; - // return FPDecimal { num: r, sign: 1 }; } if v == FPDecimal::E.num { return r + FPDecimal::ONE; - // return FPDecimal { - // num: FPDecimal::ONE.num + r, - // sign: 1, - // }; } let frac_1_5_fpdec = FPDecimal { @@ -443,7 +449,6 @@ impl FPDecimal { } pub fn log(&self, base: FPDecimal) -> FPDecimal { - // Base :self assert!(base > FPDecimal::ZERO); if *self == FPDecimal::ONE { return FPDecimal::ZERO;