Skip to content

Commit

Permalink
simplify the handling of relative tolerances
Browse files Browse the repository at this point in the history
  • Loading branch information
KRM7 committed Aug 17, 2023
1 parent b301afb commit c5a2260
Show file tree
Hide file tree
Showing 7 changed files with 47 additions and 46 deletions.
4 changes: 0 additions & 4 deletions src/utility/math.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -14,10 +14,6 @@

namespace gapp::math
{
constinit std::atomic<double> Tolerances::absolute_tolerance = 1E-12;
constinit std::atomic<unsigned> Tolerances::relative_tolerance_epsilons = 10;


bool paretoCompareLess(std::span<const double> lhs, std::span<const double> rhs) noexcept
{
GAPP_ASSERT(lhs.size() == rhs.size());
Expand Down
47 changes: 26 additions & 21 deletions src/utility/math.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -30,13 +30,13 @@ namespace gapp::math
template<std::floating_point T = double>
static T abs() noexcept { return T(absolute_tolerance.load(std::memory_order_acquire)); }

/** @returns The current relative tolerance used for floating-point comparisons. */
/** @returns The current relative tolerance used for floating-point comparisons around @p at. */
template<std::floating_point T = double>
static T eps() noexcept { return relative_tolerance_epsilons.load(std::memory_order_acquire) * std::numeric_limits<T>::epsilon(); }
static T rel(T at) noexcept { return relative_tolerance.load(std::memory_order_acquire) * at; }

private:
GAPP_API static std::atomic<double> absolute_tolerance;
GAPP_API static std::atomic<unsigned> relative_tolerance_epsilons;
GAPP_API inline constinit static std::atomic<double> absolute_tolerance = 1E-12;
GAPP_API inline constinit static std::atomic<double> relative_tolerance = 10 * std::numeric_limits<double>::epsilon();

friend class ScopedTolerances;
};
Expand All @@ -57,19 +57,22 @@ namespace gapp::math
* Create an instance of the class, setting new values for the tolerances
* used for floating-point comparisons.
*
* @param num_epsilons The number of epsilons to use as the relative tolerance in the comparisons.
* @param abs The absolute tolerance value used for the comparisons.
* @param abs The absolute tolerance value that will be used for the comparisons. Can't be negative.
* @param rel The relative tolerance value around 1.0 that will be used for the comparisons. Can't be negative.
*/
ScopedTolerances(unsigned num_epsilons, double abs) noexcept :
old_abs_tol(Tolerances::absolute_tolerance.exchange(abs, std::memory_order_acq_rel)),
old_eps_tol(Tolerances::relative_tolerance_epsilons.exchange(num_epsilons, std::memory_order_acq_rel))
{}
ScopedTolerances(double abs, double rel) noexcept :
old_absolute_tolerance(Tolerances::absolute_tolerance.exchange(abs, std::memory_order_acq_rel)),
old_relative_tolerance(Tolerances::relative_tolerance.exchange(rel, std::memory_order_acq_rel))
{
GAPP_ASSERT(abs >= 0.0, "The absolute tolerance value can't be negative.");
GAPP_ASSERT(rel >= 0.0, "The relative tolerance value can't be negative.");
}

/** Reset the tolerances to their previous values. */
~ScopedTolerances() noexcept
{
Tolerances::absolute_tolerance.store(old_abs_tol, std::memory_order_release);
Tolerances::relative_tolerance_epsilons.store(old_eps_tol, std::memory_order_release);
Tolerances::absolute_tolerance.store(old_absolute_tolerance, std::memory_order_release);
Tolerances::relative_tolerance.store(old_relative_tolerance, std::memory_order_release);
}

ScopedTolerances(const ScopedTolerances&) = delete;
Expand All @@ -78,14 +81,17 @@ namespace gapp::math
ScopedTolerances& operator=(ScopedTolerances&&) = delete;

private:
double old_abs_tol;
unsigned old_eps_tol;
double old_absolute_tolerance;
double old_relative_tolerance;
};


template<typename T>
inline constexpr T inf = std::numeric_limits<T>::infinity();

template<typename T>
inline constexpr T eps = std::numeric_limits<T>::epsilon();

template<typename T>
inline constexpr T small = std::numeric_limits<T>::min();

Expand Down Expand Up @@ -182,9 +188,8 @@ namespace gapp::math
GAPP_ASSERT(!std::isnan(lhs) && !std::isnan(rhs));

const T diff = lhs - rhs;
const T scale = std::max(std::abs(lhs), std::abs(rhs));
const T rel_tol = std::min(scale, std::numeric_limits<T>::max()) * Tolerances::eps<T>();
const T tol = std::max(rel_tol, Tolerances::abs<T>());
const T scale = std::min(std::max(std::abs(lhs), std::abs(rhs)), std::numeric_limits<T>::max());
const T tol = std::max(Tolerances::rel<T>(scale), Tolerances::abs<T>());

if (diff > tol) return 1; // lhs < rhs
if (diff < -tol) return -1; // lhs > rhs
Expand All @@ -198,7 +203,7 @@ namespace gapp::math

if (scale == inf<T>) return lhs == rhs; // for infinities

return std::abs(lhs - rhs) <= std::max(scale * Tolerances::eps<T>(), Tolerances::abs<T>());
return std::abs(lhs - rhs) <= std::max(Tolerances::rel<T>(scale), Tolerances::abs<T>());
}

template<std::floating_point T>
Expand All @@ -208,7 +213,7 @@ namespace gapp::math

if (scale == inf<T>) return lhs < rhs; // for infinities

return (rhs - lhs) > std::max(scale * Tolerances::eps<T>(), Tolerances::abs<T>());
return (rhs - lhs) > std::max(Tolerances::rel<T>(scale), Tolerances::abs<T>());
}

template<std::floating_point T>
Expand All @@ -218,7 +223,7 @@ namespace gapp::math

if (scale == inf<T>) return lhs < rhs; // for infinities

return (rhs - lhs) > std::max(scale * Tolerances::eps<T>(), Tolerances::abs<T>());
return (rhs - lhs) > std::max(Tolerances::rel<T>(scale), Tolerances::abs<T>());
}

template<std::floating_point T>
Expand All @@ -228,7 +233,7 @@ namespace gapp::math

if (scale == inf<T>) return lhs > rhs; // for infinities

return (lhs - rhs) > std::max(scale * Tolerances::eps<T>(), Tolerances::abs<T>());
return (lhs - rhs) > std::max(Tolerances::rel<T>(scale), Tolerances::abs<T>());
}

template<std::floating_point T>
Expand Down
28 changes: 14 additions & 14 deletions test/unit/math.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -22,9 +22,9 @@ TEST_CASE("fp_compare", "[math]")

SECTION("is_equal")
{
auto [rel, abs] = GENERATE(std::pair{ 0, 0.0 }, std::pair{ 10, 1E-12 });
auto [abs, rel] = GENERATE(std::pair{ 0.0, 0.0 }, std::pair{ 1E-12, 10 * eps<double> });

ScopedTolerances _(rel, abs);
ScopedTolerances _(abs, rel);
INFO("Relative tolerance eps: " << rel << ", absolute tolerance: " << abs);

REQUIRE(floatIsEqual(0.0, 0.0));
Expand Down Expand Up @@ -69,7 +69,7 @@ TEST_CASE("fp_compare", "[math]")

SECTION("approx is_equal")
{
ScopedTolerances _(10, 1E-12);
ScopedTolerances _(1E-12, 10 * eps<double>);

REQUIRE(floatIsEqual(0.0, 1E-13));
REQUIRE(!floatIsEqual(0.0, 1E-11));
Expand All @@ -80,9 +80,9 @@ TEST_CASE("fp_compare", "[math]")

SECTION("is_less")
{
auto [rel, abs] = GENERATE(std::pair{ 0, 0.0 }, std::pair{ 10, 1E-12 });
auto [abs, rel] = GENERATE(std::pair{ 0.0, 0.0 }, std::pair{ 1E-12, 10 * eps<double> });

ScopedTolerances _(rel, abs);
ScopedTolerances _(abs, rel);
INFO("Relative tolerance eps: " << rel << ", absolute tolerance: " << abs);

REQUIRE(!floatIsLess(0.0, 0.0));
Expand Down Expand Up @@ -134,7 +134,7 @@ TEST_CASE("fp_compare", "[math]")

SECTION("approx is_less")
{
ScopedTolerances _(10, 1E-12);
ScopedTolerances _(1E-12, 10 * eps<double>);

REQUIRE(!floatIsLess(0.0, 1E-13));
REQUIRE(floatIsLess(0.0, 1E-11));
Expand All @@ -145,9 +145,9 @@ TEST_CASE("fp_compare", "[math]")

SECTION("three-way comparison")
{
auto [rel, abs] = GENERATE(std::pair{ 0, 0.0 }, std::pair{ 10, 1E-12 });
auto [abs, rel] = GENERATE(std::pair{ 0.0, 0.0 }, std::pair{ 1E-12, 10 * eps<double> });

ScopedTolerances _(rel, abs);
ScopedTolerances _(abs, rel);
INFO("Relative tolerance eps: " << rel << ", absolute tolerance: " << abs);

REQUIRE(floatCompare(0.0, 0.0) == 0);
Expand Down Expand Up @@ -186,9 +186,9 @@ TEST_CASE("fp_compare", "[math]")

SECTION("is_less_not_greater")
{
auto [rel, abs] = GENERATE(std::pair{ 0, 0.0 }, std::pair{ 10, 1E-12 });
auto [abs, rel] = GENERATE(std::pair{ 0.0, 0.0 }, std::pair{ 1E-12, 10 * eps<double> });

ScopedTolerances _(rel, abs);
ScopedTolerances _(abs, rel);
INFO("Relative tolerance eps: " << rel << ", absolute tolerance: " << abs);

REQUIRE(!floatIsLessAssumeNotGreater(0.0, 0.0));
Expand Down Expand Up @@ -230,9 +230,9 @@ TEST_CASE("fp_compare", "[math]")

TEST_CASE("pareto_compare_less", "[math]")
{
auto [rel, abs] = GENERATE(std::pair{ 0, 0.0 }, std::pair{ 10, 1E-12 });
auto [abs, rel] = GENERATE(std::pair{ 0.0, 0.0 }, std::pair{ 1E-12, 10 * eps<double> });

ScopedTolerances _(rel, abs);
ScopedTolerances _(abs, rel);
INFO("Relative tolerance eps: " << rel << ", absolute tolerance: " << abs);

const std::vector vec = { 3.0, 2.0, 1.0 };
Expand Down Expand Up @@ -274,9 +274,9 @@ TEST_CASE("pareto_compare_less", "[math]")

TEST_CASE("pareto_compare_three_way", "[math]")
{
auto [rel, abs] = GENERATE(std::pair{ 0, 0.0 }, std::pair{ 10, 1E-12 });
auto [abs, rel] = GENERATE(std::pair{ 0.0, 0.0 }, std::pair{ 1E-12, 10 * eps<double> });

ScopedTolerances _(rel, abs);
ScopedTolerances _(abs, rel);
INFO("Relative tolerance eps: " << rel << ", absolute tolerance: " << abs);

const std::vector vec = { 3.0, 2.0, 1.0 };
Expand Down
2 changes: 1 addition & 1 deletion test/unit/nd_sort.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -76,7 +76,7 @@ TEMPLATE_TEST_CASE_SIG("nd_sort", "[pareto_front]", ((auto F), F), fastNonDomina
REQUIRE(std::adjacent_find(pareto_fronts.begin(), pareto_fronts.end(), [](auto lhs, auto rhs) { return (rhs.rank - lhs.rank) > 1; }) == pareto_fronts.end());


math::ScopedTolerances _(0, 0.11);
math::ScopedTolerances _(0.11, 0.0);

ParetoFronts pareto_fronts_approx = F(fmat.begin(), fmat.end());
expected_fronts[18].rank = 3;
Expand Down
8 changes: 4 additions & 4 deletions test/unit/pareto_front.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -27,7 +27,7 @@ TEST_CASE("find_pareto_front_1D", "[pareto_front]")

SECTION("multiple optimum")
{
ScopedTolerances _(0, 0.0);
ScopedTolerances _(0.0, 0.0);

auto optimal_indices = findParetoFront1D(fmat);

Expand All @@ -36,7 +36,7 @@ TEST_CASE("find_pareto_front_1D", "[pareto_front]")

SECTION("multiple optimum approx")
{
ScopedTolerances _(0, 0.1);
ScopedTolerances _(0.1, 0.0);

auto optimal_indices = findParetoFront1D(fmat);

Expand Down Expand Up @@ -75,7 +75,7 @@ TEMPLATE_TEST_CASE_SIG("find_pareto_front_nd", "[pareto_front]", ((auto F), F),

SECTION("multiple optimum")
{
ScopedTolerances _(0, 0.0);
ScopedTolerances _(0.0, 0.0);

auto optimal_indices = F(fmat);

Expand All @@ -84,7 +84,7 @@ TEMPLATE_TEST_CASE_SIG("find_pareto_front_nd", "[pareto_front]", ((auto F), F),

SECTION("multiple optimum approx")
{
ScopedTolerances _(0, 0.1);
ScopedTolerances _(0.1, 0.0);

auto optimal_indices = F(fmat);

Expand Down
2 changes: 1 addition & 1 deletion test/unit/pareto_sets.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -27,7 +27,7 @@ static constexpr bool fcomp(const Candidate<T>& lhs, const Candidate<T>& rhs) {

TEST_CASE("merge_pareto_sets", "[pareto_front]")
{
ScopedTolerances _(0, 0.0);
ScopedTolerances _(0.0, 0.0);

FitnessMatrix front1 = {
{ 10.0, -1.0 }, //
Expand Down
2 changes: 1 addition & 1 deletion test/unit/replacement.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -41,7 +41,7 @@ static const FitnessMatrix fitness_mat = {

TEST_CASE("replacement_best", "[replacement][single-objective]")
{
math::ScopedTolerances _(0, 0.0);
math::ScopedTolerances _(0.0, 0.0);

std::unique_ptr<Replacement> replacement = std::make_unique<KeepBest>();
const auto indices = replacement->nextPopulationImpl(context, fitness_mat.begin(), fitness_mat.begin() + POPSIZE, fitness_mat.end());
Expand Down

0 comments on commit c5a2260

Please sign in to comment.