Skip to content

Commit

Permalink
English.
Browse files Browse the repository at this point in the history
  • Loading branch information
@pleroy
pleroy committed Jul 19, 2024
1 parent bbbea9a commit 679ec9a
Show file tree
Hide file tree
Showing 3 changed files with 59 additions and 51 deletions.
16 changes: 8 additions & 8 deletions functions/accurate_table_generator.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -33,13 +33,13 @@ std::vector<cpp_rational> GalExhaustiveMultisearch(

// Searches in an interval of radius |T / N| centered on |starting_argument|.
// The |polynomials| must be the degree-2 Taylor approximations of the
// |functions| and the |rests| must be upper bounds on the rest of the Taylor
// series. The argument and function values must be within [1/2, 1[.
// |functions| and the |remainders| must be upper bounds on the remainder of the
// Taylor series. The argument and function values must be within [1/2, 1[.
template<std::int64_t zeroes>
absl::StatusOr<cpp_rational> StehléZimmermannSimultaneousSearch(
std::array<AccurateFunction, 2> const& functions,
std::array<AccuratePolynomial<cpp_rational, 2>, 2> const& polynomials,
std::array<AccurateFunction, 2> const& rests,
std::array<AccurateFunction, 2> const& remainders,
cpp_rational const& starting_argument,
std::int64_t N,
std::int64_t T);
Expand All @@ -51,19 +51,19 @@ template<std::int64_t zeroes>
absl::StatusOr<cpp_rational> StehléZimmermannSimultaneousFullSearch(
std::array<AccurateFunction, 2> const& functions,
std::array<AccuratePolynomial<cpp_rational, 2>, 2> const& polynomials,
std::array<AccurateFunction, 2> const& rests,
std::array<AccurateFunction, 2> const& remainders,
cpp_rational const& starting_argument);

// Same as above, but performs searches in parallel using the corresponding
// |polynomials|, |rests|, and |starting_arguments|. Returns the results in the
// same order as the parameters.
// |polynomials|, |remainders|, and |starting_arguments|. Returns the results
// in the same order as the parameters.
template<std::int64_t zeroes>
std::vector<absl::StatusOr<cpp_rational>>
StehléZimmermannSimultaneousMultisearch(
std::array<AccurateFunction, 2> const& functions,
std::vector<std::array<AccuratePolynomial<cpp_rational, 2>, 2>> const&
polynomials,
std::vector<std::array<AccurateFunction, 2>> const& rests,
std::vector<std::array<AccurateFunction, 2>> const& remainders,
std::vector<cpp_rational> const& starting_arguments);

// Same as above, but instead of accumulating all the results and returning them
Expand All @@ -74,7 +74,7 @@ void StehléZimmermannSimultaneousStreamingMultisearch(
std::array<AccurateFunction, 2> const& functions,
std::vector<std::array<AccuratePolynomial<cpp_rational, 2>, 2>> const&
polynomials,
std::vector<std::array<AccurateFunction, 2>> const& rests,
std::vector<std::array<AccurateFunction, 2>> const& remainders,
std::vector<cpp_rational> const& starting_arguments,
std::function<void(/*index=*/std::int64_t,
absl::StatusOr<cpp_rational>)> const& callback);
Expand Down
38 changes: 19 additions & 19 deletions functions/accurate_table_generator_body.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -189,7 +189,7 @@ template<std::int64_t zeroes>
absl::StatusOr<cpp_rational> StehléZimmermannSimultaneousSearch(
std::array<AccurateFunction, 2> const& functions,
std::array<AccuratePolynomial<cpp_rational, 2>, 2> const& polynomials,
std::array<AccurateFunction, 2> const& rests,
std::array<AccurateFunction, 2> const& remainders,
cpp_rational const& starting_argument,
std::int64_t const N,
std::int64_t const T) {
Expand Down Expand Up @@ -225,14 +225,14 @@ absl::StatusOr<cpp_rational> StehléZimmermannSimultaneousSearch(
P[i] = N * Compose(polynomials[i], shift_and_rescale);
}

// Step 2: compute ε. We use the rests provided by the clients. Note that we
// could save on the number of evaluations by providing both bounds to a
// single call.
// Step 2: compute ε. We use the remainders provided by the clients. Note
// that we could save on the number of evaluations by providing both bounds to
// a single call.
cpp_bin_float_50 ε = 0;
for (std::int64_t i = 0; i < rests.size(); ++i) {
for (std::int64_t i = 0; i < remainders.size(); ++i) {
auto const T_over_N = cpp_rational(T, N);
ε = std::max(ε, abs(N * rests[i](starting_argument - T_over_N)));
ε = std::max(ε, abs(N * rests[i](starting_argument + T_over_N)));
ε = std::max(ε, abs(N * remainders[i](starting_argument - T_over_N)));
ε = std::max(ε, abs(N * remainders[i](starting_argument + T_over_N)));
}
VLOG(2) << "ε = " << ε;

Expand Down Expand Up @@ -371,7 +371,7 @@ template<std::int64_t zeroes>
absl::StatusOr<cpp_rational> StehléZimmermannSimultaneousFullSearch(
std::array<AccurateFunction, 2> const& functions,
std::array<AccuratePolynomial<cpp_rational, 2>, 2> const& polynomials,
std::array<AccurateFunction, 2> const& rests,
std::array<AccurateFunction, 2> const& remainders,
cpp_rational const& starting_argument) {
constexpr std::int64_t M = 1LL << zeroes;
constexpr std::int64_t N = 1LL << std::numeric_limits<double>::digits;
Expand All @@ -392,7 +392,7 @@ absl::StatusOr<cpp_rational> StehléZimmermannSimultaneousFullSearch(

std::array<double, 2> function_scales;
std::array<AccurateFunction, 2> scaled_functions;
std::array<AccurateFunction, 2> scaled_rests;
std::array<AccurateFunction, 2> scaled_remainders;
for (std::int64_t i = 0; i < scaled_functions.size(); ++i) {
std::int64_t function_exponent;
auto const function_mantissa =
Expand All @@ -405,11 +405,11 @@ absl::StatusOr<cpp_rational> StehléZimmermannSimultaneousFullSearch(
i](cpp_rational const& argument) {
return function_scale * functions[i](argument / argument_scale);
};
scaled_rests[i] = [argument_scale,
scaled_remainders[i] = [argument_scale,
function_scale = function_scales[i],
&rests,
&remainders,
i](cpp_rational const& argument) {
return function_scale * rests[i](argument / argument_scale);
return function_scale * remainders[i](argument / argument_scale);
};
}

Expand Down Expand Up @@ -442,7 +442,7 @@ absl::StatusOr<cpp_rational> StehléZimmermannSimultaneousFullSearch(
auto const status_or_solution =
StehléZimmermannSimultaneousSearch<zeroes>(scaled_functions,
scaled_polynomials,
scaled_rests,
scaled_remainders,
high_interval.midpoint(),
N,
T);
Expand Down Expand Up @@ -479,7 +479,7 @@ absl::StatusOr<cpp_rational> StehléZimmermannSimultaneousFullSearch(
auto const status_or_solution =
StehléZimmermannSimultaneousSearch<zeroes>(scaled_functions,
scaled_polynomials,
scaled_rests,
scaled_remainders,
low_interval.midpoint(),
N,
T);
Expand Down Expand Up @@ -526,14 +526,14 @@ StehléZimmermannSimultaneousMultisearch(
std::array<AccurateFunction, 2> const& functions,
std::vector<std::array<AccuratePolynomial<cpp_rational, 2>, 2>> const&
polynomials,
std::vector<std::array<AccurateFunction, 2>> const& rests,
std::vector<std::array<AccurateFunction, 2>> const& remainders,
std::vector<cpp_rational> const& starting_arguments) {
std::vector<absl::StatusOr<cpp_rational>> result;
result.resize(starting_arguments.size());
StehléZimmermannSimultaneousStreamingMultisearch<zeroes>(
functions,
polynomials,
rests,
remainders,
starting_arguments,
[&result](std::int64_t const index,
absl::StatusOr<cpp_rational> status_or_final_argument) {
Expand All @@ -547,7 +547,7 @@ void StehléZimmermannSimultaneousStreamingMultisearch(
std::array<AccurateFunction, 2> const& functions,
std::vector<std::array<AccuratePolynomial<cpp_rational, 2>, 2>> const&
polynomials,
std::vector<std::array<AccurateFunction, 2>> const& rests,
std::vector<std::array<AccurateFunction, 2>> const& remainders,
std::vector<cpp_rational> const& starting_arguments,
std::function<void(/*index=*/std::int64_t,
absl::StatusOr<cpp_rational>)> const& callback) {
Expand All @@ -559,13 +559,13 @@ void StehléZimmermannSimultaneousStreamingMultisearch(
&callback,
&functions,
&polynomials,
&rests,
&remainders,
&starting_arguments]() {
auto const& starting_argument = starting_arguments[i];
LOG(INFO) << "Starting search around " << starting_argument;
auto status_or_final_argument =
StehléZimmermannSimultaneousFullSearch<zeroes>(
functions, polynomials[i], rests[i], starting_argument);
functions, polynomials[i], remainders[i], starting_argument);
if (status_or_final_argument.ok()) {
LOG(INFO) << "Finished search around " << starting_argument
<< ", found " << status_or_final_argument.value();
Expand Down
56 changes: 32 additions & 24 deletions functions/accurate_table_generator_test.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -143,14 +143,16 @@ TEST_F(AccurateTableGeneratorTest, StehléZimmermannSinCos15) {
-cpp_rational(Cos(u₀ / 4) / 32)},
u₀);

// Use the Lagrange form of the rest. u may be above or below u₀. Note that
// we use the fact that the functions are monotonic.
auto const rest_sin_taylor2 = [u₀ = cpp_rational(u₀)](cpp_rational const& u) {
// Use the Lagrange form of the remainder. u may be above or below u₀. Note
// that we use the fact that the functions are monotonic.
auto const remainder_sin_taylor2 =
[u₀ = cpp_rational(u₀)](cpp_rational const& u) {
auto const Δu = u - u₀;
auto const Δu³ = static_cast<cpp_bin_float_50>(Δu * Δu * Δu);
return Δu³ * (-Cos(std::min(u, u₀) / 4) / 16) / Factorial(3);
};
auto const rest_cos_taylor2 = [u₀ = cpp_rational(u₀)](cpp_rational const& u) {
auto const remainder_cos_taylor2 =
[u₀ = cpp_rational(u₀)](cpp_rational const& u) {
auto const Δu = u - u₀;
auto const Δu³ = static_cast<cpp_bin_float_50>(Δu * Δu * Δu);
return Δu³ * (Sin(std::max(u, u₀) / 4) / 64) / Factorial(3);
Expand All @@ -159,7 +161,7 @@ TEST_F(AccurateTableGeneratorTest, StehléZimmermannSinCos15) {
auto const u = StehléZimmermannSimultaneousSearch<15>(
{sin, cos},
{sin_taylor2, cos_taylor2},
{rest_sin_taylor2, rest_cos_taylor2},
{remainder_sin_taylor2, remainder_cos_taylor2},
u₀,
/*N=*/1ll << 53,
/*T=*/1ll << 21);
Expand Down Expand Up @@ -204,12 +206,14 @@ TEST_F(AccurateTableGeneratorTest, StehléZimmermannFullSinCos5NoScaling) {
-cpp_rational(Cos(u₀ / 4) / 32)},
u₀);

auto const rest_sin_taylor2 = [u₀ = cpp_rational(u₀)](cpp_rational const& u) {
auto const remainder_sin_taylor2 =
[u₀ = cpp_rational(u₀)](cpp_rational const& u) {
auto const Δu = u - u₀;
auto const Δu³ = static_cast<cpp_bin_float_50>(Δu * Δu * Δu);
return Δu³ * (-Cos(std::min(u, u₀) / 4) / 16) / Factorial(3);
};
auto const rest_cos_taylor2 = [u₀ = cpp_rational(u₀)](cpp_rational const& u) {
auto const remainder_cos_taylor2 =
[u₀ = cpp_rational(u₀)](cpp_rational const& u) {
auto const Δu = u - u₀;
auto const Δu³ = static_cast<cpp_bin_float_50>(Δu * Δu * Δu);
return Δu³ * (Sin(std::max(u, u₀) / 4) / 64) / Factorial(3);
Expand All @@ -218,7 +222,7 @@ TEST_F(AccurateTableGeneratorTest, StehléZimmermannFullSinCos5NoScaling) {
auto const u = StehléZimmermannSimultaneousFullSearch<5>(
{sin, cos},
{sin_taylor2, cos_taylor2},
{rest_sin_taylor2, rest_cos_taylor2},
{remainder_sin_taylor2, remainder_cos_taylor2},
u₀);
EXPECT_THAT(u,
IsOkAndHolds(cpp_rational(4785074604080979, 9007199254740992)));
Expand Down Expand Up @@ -261,12 +265,14 @@ TEST_F(AccurateTableGeneratorTest, StehléZimmermannFullSinCos15NoScaling) {
-cpp_rational(Cos(u₀ / 4) / 32)},
u₀);

auto const rest_sin_taylor2 = [u₀ = cpp_rational(u₀)](cpp_rational const& u) {
auto const remainder_sin_taylor2 =
[u₀ = cpp_rational(u₀)](cpp_rational const& u) {
auto const Δu = u - u₀;
auto const Δu³ = static_cast<cpp_bin_float_50>(Δu * Δu * Δu);
return Δu³ * (-Cos(std::min(u, u₀) / 4) / 16) / Factorial(3);
};
auto const rest_cos_taylor2 = [u₀ = cpp_rational(u₀)](cpp_rational const& u) {
auto const remainder_cos_taylor2 =
[u₀ = cpp_rational(u₀)](cpp_rational const& u) {
auto const Δu = u - u₀;
auto const Δu³ = static_cast<cpp_bin_float_50>(Δu * Δu * Δu);
return Δu³ * (Sin(std::max(u, u₀) / 4) / 64) / Factorial(3);
Expand All @@ -275,7 +281,7 @@ TEST_F(AccurateTableGeneratorTest, StehléZimmermannFullSinCos15NoScaling) {
auto const u = StehléZimmermannSimultaneousFullSearch<15>(
{sin, cos},
{sin_taylor2, cos_taylor2},
{rest_sin_taylor2, rest_cos_taylor2},
{remainder_sin_taylor2, remainder_cos_taylor2},
u₀);
EXPECT_THAT(u,
IsOkAndHolds(cpp_rational(4785074575333183, 9007199254740992)));
Expand Down Expand Up @@ -310,12 +316,14 @@ TEST_F(AccurateTableGeneratorTest, StehléZimmermannFullSinCos15WithScaling) {
-cpp_rational(Cos(x₀) / 2)},
x₀);

auto const rest_sin_taylor2 = [x₀ = cpp_rational(x₀)](cpp_rational const& x) {
auto const remainder_sin_taylor2 =
[x₀ = cpp_rational(x₀)](cpp_rational const& x) {
auto const Δx = x - x₀;
auto const Δx³ = static_cast<cpp_bin_float_50>(Δx * Δx * Δx);
return -Δx³ * -Cos(std::min(x, x₀)) / Factorial(3);
};
auto const rest_cos_taylor2 = [x₀ = cpp_rational(x₀)](cpp_rational const& x) {
auto const remainder_cos_taylor2 =
[x₀ = cpp_rational(x₀)](cpp_rational const& x) {
auto const Δx = x - x₀;
auto const Δx³ = static_cast<cpp_bin_float_50>(Δx * Δx * Δx);
return Δx³ * Sin(std::max(x, x₀)) / Factorial(3);
Expand All @@ -324,7 +332,7 @@ TEST_F(AccurateTableGeneratorTest, StehléZimmermannFullSinCos15WithScaling) {
auto const x = StehléZimmermannSimultaneousFullSearch<15>(
{Sin, Cos},
{sin_taylor2, cos_taylor2},
{rest_sin_taylor2, rest_cos_taylor2},
{remainder_sin_taylor2, remainder_cos_taylor2},
x₀);
EXPECT_THAT(x,
IsOkAndHolds(cpp_rational(4785074575333183, 36028797018963968)));
Expand Down Expand Up @@ -353,7 +361,7 @@ TEST_F(AccurateTableGeneratorTest, StehléZimmermannMultisearchSinCos15) {
static constexpr std::int64_t index_end = 100;
std::vector<cpp_rational> starting_arguments;
std::vector<std::array<AccuratePolynomial<cpp_rational, 2>, 2>> polynomials;
std::vector<std::array<AccurateFunction, 2>> rests;
std::vector<std::array<AccurateFunction, 2>> remainders;
for (std::int64_t i = index_begin; i < index_end; ++i) {
auto const x₀ = i / 128.0;
AccuratePolynomial<cpp_rational, 2> const sin_taylor2(
Expand All @@ -367,13 +375,13 @@ TEST_F(AccurateTableGeneratorTest, StehléZimmermannMultisearchSinCos15) {
-cpp_rational(Cos(x₀) / 2)},
x₀);

auto const rest_sin_taylor2 =
auto const remainder_sin_taylor2 =
[x₀ = cpp_rational(x₀)](cpp_rational const& x) {
auto const Δx = x - x₀;
auto const Δx³ = static_cast<cpp_bin_float_50>(Δx * Δx * Δx);
return -Δx³ * Cos(std::min(x, x₀)) / Factorial(3);
};
auto const rest_cos_taylor2 =
auto const remainder_cos_taylor2 =
[x₀ = cpp_rational(x₀)](cpp_rational const& x) {
auto const Δx = x - x₀;
auto const Δx³ = static_cast<cpp_bin_float_50>(Δx * Δx * Δx);
Expand All @@ -382,10 +390,10 @@ TEST_F(AccurateTableGeneratorTest, StehléZimmermannMultisearchSinCos15) {

starting_arguments.push_back(x₀);
polynomials.push_back({sin_taylor2, cos_taylor2});
rests.push_back({rest_sin_taylor2, rest_cos_taylor2});
remainders.push_back({remainder_sin_taylor2, remainder_cos_taylor2});
}
auto const xs = StehléZimmermannSimultaneousMultisearch<15>(
{Sin, Cos}, polynomials, rests, starting_arguments);
{Sin, Cos}, polynomials, remainders, starting_arguments);
EXPECT_THAT(xs, SizeIs(index_end - index_begin));
for (std::int64_t i = 0; i < xs.size(); ++i) {
CHECK_OK(xs[i].status());
Expand Down Expand Up @@ -430,7 +438,7 @@ TEST_F(AccurateTableGeneratorTest, DISABLED_SECULAR_SinCos18) {

std::vector<cpp_rational> starting_arguments;
std::vector<std::array<AccuratePolynomial<cpp_rational, 2>, 2>> polynomials;
std::vector<std::array<AccurateFunction, 2>> rests;
std::vector<std::array<AccurateFunction, 2>> remainders;
for (std::int64_t i = std::floor(lower_bound / h_over_2);
i <= std::ceil(upper_bound / h_over_2);
++i) {
Expand All @@ -449,13 +457,13 @@ TEST_F(AccurateTableGeneratorTest, DISABLED_SECULAR_SinCos18) {
-cpp_rational(Cos(x₀) / 2)},
x₀);

auto const rest_sin_taylor2 =
auto const remainder_sin_taylor2 =
[x₀ = cpp_rational(x₀)](cpp_rational const& x) {
auto const Δx = x - x₀;
auto const Δx³ = static_cast<cpp_bin_float_50>(Δx * Δx * Δx);
return -Δx³ * -Cos(std::min(x, x₀)) / Factorial(3);
};
auto const rest_cos_taylor2 =
auto const remainder_cos_taylor2 =
[x₀ = cpp_rational(x₀)](cpp_rational const& x) {
auto const Δx = x - x₀;
auto const Δx³ = static_cast<cpp_bin_float_50>(Δx * Δx * Δx);
Expand All @@ -464,15 +472,15 @@ TEST_F(AccurateTableGeneratorTest, DISABLED_SECULAR_SinCos18) {

starting_arguments.push_back(x₀);
polynomials.push_back({sin_taylor2, cos_taylor2});
rests.push_back({rest_sin_taylor2, rest_cos_taylor2});
remainders.push_back({remainder_sin_taylor2, remainder_cos_taylor2});
}
}
}

StehléZimmermannSimultaneousStreamingMultisearch<18>(
{Sin, Cos},
polynomials,
rests,
remainders,
starting_arguments,
[n, &logger](std::int64_t const index,
absl::StatusOr<cpp_rational> status_or_x) {
Expand Down

0 comments on commit 679ec9a

Please sign in to comment.