Skip to content

Commit

Permalink
a
Browse files Browse the repository at this point in the history
  • Loading branch information
@fintarin
fintarin committed Apr 21, 2024
1 parent d501166 commit 2a2f4c9
Show file tree
Hide file tree
Showing 4 changed files with 174 additions and 47 deletions.
9 changes: 3 additions & 6 deletions include/fintamath/numbers/IntegerFunctions.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -13,18 +13,15 @@ namespace fintamath {

using FactorToCountMap = std::unordered_map<Integer, Integer>;

// Use exponentiation by squaring with constant auxiliary memory (iterative version).
// https://en.wikipedia.org/wiki/Exponentiation_by_squaring#With_constant_auxiliary_memory.
template <std::derived_from<INumber> Lhs>
Lhs pow(const Lhs &lhs, Integer rhs) {
if (lhs == 0 && rhs == 0) {
throw UndefinedException("Undefined 0^0");
}

if (rhs < 0) {
return pow(1 / lhs, -rhs);
}

// Use exponentiation by squaring with constant auxiliary memory (iterative version).
// https://en.wikipedia.org/wiki/Exponentiation_by_squaring#With_constant_auxiliary_memory.

Lhs res(1);
Lhs sqr = lhs;

Expand Down
79 changes: 56 additions & 23 deletions src/fintamath/numbers/RealFunctions.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -2,6 +2,8 @@

#include <stdexcept>

#include <fmt/core.h>

#include <boost/math/policies/error_handling.hpp>
#include <boost/multiprecision/detail/default_ops.hpp>
#include <boost/multiprecision/mpfr.hpp>
Expand All @@ -15,6 +17,12 @@ namespace fintamath {

using namespace detail;

namespace {

bool isNegZero(const Real &rhs) {
return rhs.isZero() && rhs.sign() < 0;
}

bool isOverflow(const Real &rhs) {
static Cache<unsigned, Real::Backend> cache([](const unsigned precision) {
static const Real::Backend powBase = 10;
Expand All @@ -41,76 +49,101 @@ bool isLogUnderflow(const Real &rhs) {
return rhs != 0 && abs(rhs) < cache[Real::getPrecision()];
}

}

const Real &trigResultChecked(const Real &rhs) {
if (isUnderflow(rhs) || isOverflow(rhs)) {
throw UndefinedException("...");
if (isUnderflow(rhs)) {
throw UndefinedException("underflow");
}

if (isOverflow(rhs)) {
throw UndefinedException("overflow");
}

return rhs;
}

const Real &hyperbResultChecked(const Real &rhs) {
if (isUnderflow(rhs)) {
throw UndefinedException("...");
throw UndefinedException("underflow");
}

return rhs;
}

const Real &tgammaResultChecked(const Real &rhs) {
if (rhs == 0) {
throw UndefinedException("...");
throw UndefinedException("overflow");
}

return rhs;
}

Integer floor(const Real &rhs) {
Integer floor(const Real &rhs) try {
if (isOverflow(rhs)) {
throw UndefinedFunctionException("floor", {rhs.toString()});
throw UndefinedException("overflow");
}

const Real::Backend res = boost::multiprecision::floor(rhs.getBackend());
return res.convert_to<Integer::Backend>();
}
catch (const UndefinedException &exc) {
throw UndefinedException(fmt::format(
"Undefined floor({}) ({})",
rhs.toString(),
exc.what()));
}

Integer ceil(const Real &rhs) {
Integer ceil(const Real &rhs) try {
if (isOverflow(rhs)) {
throw UndefinedFunctionException("ceil", {rhs.toString()});
throw UndefinedException("overflow");
}

const Real::Backend res = boost::multiprecision::ceil(rhs.getBackend());
return res.convert_to<Integer::Backend>();
}
catch (const UndefinedException &exc) {
throw UndefinedException(fmt::format(
"Undefined ceil({}) ({})",
rhs.toString(),
exc.what()));
}

Real abs(const Real &rhs) {
return rhs < 0 ? -rhs : rhs;
}

Real sqrt(const Real &rhs) {
Real sqrt(const Real &rhs) try {
if (rhs < 0) {
throw UndefinedFunctionException("sqrt", {rhs.toString()});
throw UndefinedException("expected rhs >= 0");
}

return {sqrt(rhs.getBackend())};
}
catch (const UndefinedException &exc) {
throw UndefinedException(fmt::format(
"Undefined sqrt({}) ({})",
rhs.toString(),
exc.what()));
}

Real pow(const Real &lhs, const Real &rhs) {
if (abs(lhs) == 0 && abs(rhs) == 0) {
throw UndefinedBinaryOperatorException("^", lhs.toString(), rhs.toString());
Real pow(const Real &lhs, const Real &rhs) try {
if (isNegZero(lhs)) {
(void)Real(pow(-1, rhs.getBackend())); // Use (-1)^rhs to validate (-0)^rhs
}

try {
if (lhs == Real("-0")) {
// Use (-1)^rhs to validate (-0)^rhs
Real(pow(-1, rhs.getBackend()));
}

return {pow(lhs.getBackend(), rhs.getBackend())};
}
catch (const UndefinedException &) {
throw UndefinedBinaryOperatorException("^", lhs.toString(), rhs.toString());
if (isNegZero(rhs)) {
(void)Real(pow(lhs.getBackend(), -1)); // Use lhs^(-1) to validate lhs^(-0)
}

return {pow(lhs.getBackend(), rhs.getBackend())};
}
catch (const UndefinedException &exc) {
throw UndefinedException(fmt::format(
"Undefined pow({}, {}) ({})",
lhs.toString(),
rhs.toString(),
exc.what()));
}

Real exp(const Real &rhs) {
Expand Down
25 changes: 8 additions & 17 deletions tests/src/numbers/IntegerFunctionsTests.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -124,6 +124,8 @@ TEST(IntegerFunctionsTests, sqrtWithRemainderTest) {
}

TEST(IntegerFunctionsTests, powTest) {
EXPECT_EQ(pow(Integer(0), Integer(0)),
1);
EXPECT_EQ(pow(Integer(5), Integer(2)),
25);
EXPECT_EQ(pow(Integer(-5), Integer(5)),
Expand All @@ -137,6 +139,8 @@ TEST(IntegerFunctionsTests, powTest) {
EXPECT_EQ(pow(Integer(6789), Integer(-4)),
0);

EXPECT_EQ(pow(Rational(0), Integer(0)),
1);
EXPECT_EQ(pow(Rational(5, 2), Integer(2)).toString(),
"25/4");
EXPECT_EQ(pow(Rational(-5, 2), Integer(5)).toString(),
Expand All @@ -148,6 +152,8 @@ TEST(IntegerFunctionsTests, powTest) {
EXPECT_EQ(pow(Rational(6789), Integer(-4)).toString(),
"1/2124336126051441");

EXPECT_EQ(pow(Real(0), Integer(0)),
1);
EXPECT_EQ(pow(Real(5), Integer(2)).toString(),
"25.0");
EXPECT_EQ(pow(Real(-5), Integer(5)).toString(),
Expand All @@ -159,6 +165,8 @@ TEST(IntegerFunctionsTests, powTest) {
EXPECT_EQ(pow(Real(6789), Integer(-4)).toString(),
"4.7073529830645308412*10^-16");

EXPECT_EQ(pow(Complex(0), Integer(0)),
1);
EXPECT_EQ(pow(Complex(5, 2), Integer(2)).toString(),
"21 + 20 I");
EXPECT_EQ(pow(Complex(Rational(1, 2), Rational(2, 3)), Integer(5)).toString(),
Expand All @@ -169,23 +177,6 @@ TEST(IntegerFunctionsTests, powTest) {
"4832537992678386348337867205980822373798257851094432575433371642956047093945525418683507600681355491343220122262032882136893556623485326120689398001084489/1953125000000000000000000000000000000000000000000");
EXPECT_EQ(pow(Complex(6789, 11), Integer(-4)).toString(),
"531075666086959/1128212841481282934557153710724 - 860496245400/282053210370320733639288427681 I");

EXPECT_THAT(
[] { pow(Integer(0), Integer(0)); },
testing::ThrowsMessage<UndefinedException>(
testing::StrEq(R"(Undefined 0^0)")));
EXPECT_THAT(
[] { pow(Rational(0), Integer(0)); },
testing::ThrowsMessage<UndefinedException>(
testing::StrEq(R"(Undefined 0^0)")));
EXPECT_THAT(
[] { pow(Real(0), Integer(0)); },
testing::ThrowsMessage<UndefinedException>(
testing::StrEq(R"(Undefined 0^0)")));
EXPECT_THAT(
[] { pow(Complex(0), Integer(0)); },
testing::ThrowsMessage<UndefinedException>(
testing::StrEq(R"(Undefined 0^0)")));
}

TEST(IntegerFunctionsTests, factorialTest) {
Expand Down
Loading

0 comments on commit 2a2f4c9

Please sign in to comment.