Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Example 011 #276

Merged
merged 2 commits into from
Nov 1, 2023
Merged

Example 011 #276

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
c4cb9b0
Refactor example eleven
ckormanyos Nov 1, 2023
cfd3993
Update CI script
ckormanyos Nov 1, 2023
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
10 changes: 5 additions & 5 deletions .github/workflows/wide_decimal.yml
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -51,7 +51,7 @@ jobs:
- name: test
working-directory: build
run: ctest --verbose --output-on-failure
gnumake-clang-tidy-12-native:
gnumake-clang-tidy-native:
runs-on: ubuntu-latest
defaults:
run:
Expand All @@ -77,7 +77,7 @@ jobs:
git submodule update --init libs/multiprecision
./bootstrap.sh
./b2 headers
- name: gnumake-clang-tidy-12-native
- name: gnumake-clang-tidy-native
run: |
grep BOOST_VERSION ../boost-root/boost/version.hpp
cd .tidy/make
Expand Down Expand Up @@ -451,7 +451,7 @@ jobs:
git submodule update --init libs/multiprecision
./bootstrap.bat
./b2 headers
- uses: actions/checkout@v1
- uses: actions/checkout@v3
- uses: ilammy/msvc-dev-cmd@v1
with:
toolset: 14.2
Expand Down Expand Up @@ -479,7 +479,7 @@ jobs:
git submodule update --init libs/multiprecision
./bootstrap.bat
./b2 headers
- uses: actions/checkout@v1
- uses: actions/checkout@v3
- uses: ilammy/msvc-dev-cmd@v1
with:
toolset: 14.3
Expand Down Expand Up @@ -507,7 +507,7 @@ jobs:
git submodule update --init libs/multiprecision
./bootstrap.bat
./b2 headers
- uses: actions/checkout@v1
- uses: actions/checkout@v3
- uses: ilammy/msvc-dev-cmd@v1
with:
toolset: 14.3
Expand Down
192 changes: 106 additions & 86 deletions examples/example011_trig_trapezoid_integral.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -101,7 +101,7 @@ namespace example011_trig
* x2 + coef_sin_bot_1)
* x2 + coef_sin_bot_0);

return (x * top) / (bot * 153U); // NOLINT(cppcoreguidelines-avoid-magic-numbers,readability-magic-numbers)
return (x * top) / (bot * static_cast<unsigned>(UINT8_C(153)));
}

auto cos_pade(const dec51_t& x) -> dec51_t
Expand Down Expand Up @@ -169,18 +169,18 @@ namespace example011_trig
* x2 + coef_cos_bot_1)
* x2 + coef_cos_bot_0);

return 1U + (((x2 * 4080U) * top) / bot); // NOLINT(cppcoreguidelines-avoid-magic-numbers,readability-magic-numbers)
return static_cast<unsigned>(UINT8_C(1)) + (((x2 * static_cast<unsigned>(UINT16_C(4080))) * top) / bot);
}

auto sin(const dec51_t& x) -> dec51_t // NOLINT(misc-no-recursion)
{
dec51_t s;
dec51_t s { static_cast<unsigned>(UINT8_C(0)) };

if(x < 0)
if(x < static_cast<int>(INT8_C(0)))
{
s = -sin(-x);
}
else if(x > 0)
else if(x > static_cast<int>(INT8_C(0)))
{
// Perform argument reduction and subsequent scaling of the result.

Expand All @@ -193,60 +193,65 @@ namespace example011_trig
// | 2 | -sin(r) | -cos(r) | sin(r)/cos(r) |
// | 3 | -cos(r) | sin(r) | -cos(r)/sin(r) |

const auto my_pi_half = pi<dec51_t>() / 2;
static const auto my_pi_half = pi<dec51_t>() / static_cast<unsigned>(UINT8_C(2));

const auto k = static_cast<std::uint32_t> (x / my_pi_half);
const auto n = static_cast<std::uint_fast32_t>(k % 4U);
const auto k = static_cast<unsigned>(x / my_pi_half);
const auto n = static_cast<unsigned>(k % static_cast<unsigned>(UINT8_C(4)));

dec51_t r = x - (my_pi_half * k);

const auto is_neg = (n > 1U);
const auto is_cos = ((n == 1U) || (n == 3U));
auto n_angle_identity = static_cast<unsigned>(UINT8_C(0));

auto n_angle_identity = static_cast<std::uint_fast32_t>(0U);

static const dec51_t two_tenths = dec51_t(2U) / 10U;
static const auto two_tenths = dec51_t(static_cast<unsigned>(UINT8_C(2))) / static_cast<unsigned>(UINT8_C(10));

// Reduce the argument with factors of three until it is less than 2/10.
while(r > two_tenths) // NOLINT(altera-id-dependent-backward-branch)
{
r /= 3U;
r /= static_cast<unsigned>(UINT8_C(3));

++n_angle_identity;
}

s = (is_cos ? cos_pade(r) : sin_pade(r));
switch(n)
{
case static_cast<unsigned>(UINT8_C(1)):
case static_cast<unsigned>(UINT8_C(3)):
s = cos_pade(r);
break;
case static_cast<unsigned>(UINT8_C(0)):
case static_cast<unsigned>(UINT8_C(2)):
default:
s = sin_pade(r);
break;
}

// Rescale the result with the triple angle identity for sine.
for(auto t = static_cast<std::uint_fast32_t>(0U); t < n_angle_identity; ++t)
for(auto t = static_cast<unsigned>(UINT8_C(0)); t < n_angle_identity; ++t)
{
s = (s * 3U) - (((s * s) * s) * 4U);
static_cast<void>(t);

s = (s * static_cast<unsigned>(UINT8_C(3))) - (((s * s) * s) * static_cast<unsigned>(UINT8_C(4)));
}

s = fabs(s);
if(s.isneg()) { s.negate(); }

if(is_neg)
{
s = -s;
}
}
else
{
s = 0;
const auto b_neg = (n > static_cast<unsigned>(UINT8_C(1)));

if(b_neg) { s.negate(); }
}

return s;
}

auto cos(const dec51_t& x) -> dec51_t // NOLINT(misc-no-recursion)
{
dec51_t c;
dec51_t c { static_cast<unsigned>(UINT8_C(1)) };

if(x < 0)
if(x < static_cast<int>(INT8_C(0)))
{
c = cos(-x);
}
else if(x > 0)
else if(x > static_cast<int>(INT8_C(0)))
{
// Perform argument reduction and subsequent scaling of the result.

Expand All @@ -259,46 +264,52 @@ namespace example011_trig
// | 2 | -sin(r) | -cos(r) | sin(r)/cos(r) |
// | 3 | -cos(r) | sin(r) | -cos(r)/sin(r) |

const auto my_pi_half = pi<dec51_t>() / 2;
static const auto my_pi_half = pi<dec51_t>() / static_cast<unsigned>(UINT8_C(2));

const auto k = static_cast<std::uint32_t> (x / my_pi_half);
const auto n = static_cast<std::uint_fast32_t>(k % 4U);
const auto k = static_cast<unsigned>(x / my_pi_half);
const auto n = static_cast<unsigned>(k % static_cast<unsigned>(UINT8_C(4)));

dec51_t r = x - (my_pi_half * k);

const auto is_neg = ((n == 1U) || (n == 2U));
const auto is_sin = ((n == 1U) || (n == 3U));
auto n_angle_identity = static_cast<unsigned>(UINT8_C(0));

auto n_angle_identity = static_cast<std::uint_fast32_t>(0U);

static const dec51_t two_tenths = dec51_t(2U) / 10U;
static const auto two_tenths = dec51_t(static_cast<unsigned>(UINT8_C(2))) / static_cast<unsigned>(UINT8_C(10));

// Reduce the argument with factors of three until it is less than 2/10.
while(r > two_tenths) // NOLINT(altera-id-dependent-backward-branch)
{
r /= 3U;
r /= static_cast<unsigned>(UINT8_C(3));

++n_angle_identity;
}

c = (is_sin ? sin_pade(r) : cos_pade(r));

// Rescale the result with the triple angle identity for cosine.
for(auto t = static_cast<std::uint_fast32_t>(0U); t < n_angle_identity; ++t)
switch(n)
{
c = (((c * c) * c) * 4U) - (c * 3U);
case static_cast<unsigned>(UINT8_C(1)):
case static_cast<unsigned>(UINT8_C(3)):
c = sin_pade(r);
break;
case static_cast<unsigned>(UINT8_C(0)):
case static_cast<unsigned>(UINT8_C(2)):
default:
c = cos_pade(r);
break;
}

c = fabs(c);

if(is_neg)
// Rescale the result with the triple angle identity for cosine.
for(auto t = static_cast<unsigned>(UINT8_C(0)); t < n_angle_identity; ++t)
{
c = -c;
static_cast<void>(t);

c = (((c * c) * c) * static_cast<unsigned>(UINT8_C(4))) - (c * static_cast<unsigned>(UINT8_C(3)));
}
}
else
{
c = dec51_t(1U);

if(c.isneg()) { c.negate(); }

const auto b_neg = ((n == static_cast<unsigned>(UINT8_C(1))) || (n == static_cast<unsigned>(UINT8_C(2))));

if(b_neg) { c.negate(); }
}

return c;
Expand All @@ -313,17 +324,17 @@ namespace example011_trig
{
auto n2 = static_cast<std::uint_fast32_t>(INT8_C(1));

real_value_type step = ((b - a) / 2U);
real_value_type step = ((b - a) / static_cast<unsigned>(UINT8_C(2)));

real_value_type result = (real_function(a) + real_function(b)) * step;

const auto k_max = static_cast<std::uint_fast8_t>(32U);
const auto k_max = static_cast<std::uint_fast8_t>(UINT8_C(32));

for(auto k = static_cast<std::uint_fast8_t>(0U); k < k_max; ++k)
for(auto k = static_cast<std::uint_fast8_t>(UINT8_C(0)); k < k_max; ++k)
{
real_value_type sum(0);
real_value_type sum(static_cast<unsigned>(UINT8_C(0)));

for(auto j = static_cast<std::uint_fast32_t>(0U); j < n2; ++j)
for(auto j = static_cast<std::uint_fast32_t>(UINT8_C(0)); j < n2; ++j)
{
const auto two_j_plus_one =
static_cast<std::uint_fast32_t>
Expand All @@ -336,22 +347,22 @@ namespace example011_trig

const real_value_type tmp = result;

result = (result / 2U) + (step * sum);
result = (result / static_cast<unsigned>(UINT8_C(2))) + (step * sum);

using std::fabs;

const real_value_type ratio = fabs(tmp / result);

const real_value_type delta = fabs(ratio - 1U);
const real_value_type delta = fabs(ratio - static_cast<unsigned>(UINT8_C(1)));

if((k > UINT8_C(1)) && (delta < tol))
if((k > static_cast<std::uint_fast8_t>(UINT8_C(1))) && (delta < tol))
{
break;
}

n2 *= 2U;
n2 *= static_cast<unsigned>(UINT8_C(2));

step /= 2U;
step /= static_cast<unsigned>(UINT8_C(2));
}

return result;
Expand All @@ -361,26 +372,31 @@ namespace example011_trig
auto cyl_bessel_j( std::uint_fast8_t n,
const float_type& x) -> float_type
{
const float_type epsilon = std::numeric_limits<float_type>::epsilon();

using std::cos;
using std::sin;
using std::sqrt;

using example011_trig::cos;
using example011_trig::sin;
static const auto tol = sqrt(std::numeric_limits<float_type>::epsilon());

static const auto my_pi = pi<float_type>();

const float_type integration_result =
integral
(
float_type(static_cast<unsigned>(UINT8_C(0))),
my_pi,
tol,
[&x, &n](const float_type& t) -> float_type
{
using std::cos;
using std::sin;

const float_type tol = sqrt(epsilon);
using example011_trig::cos;
using example011_trig::sin;

const float_type integration_result = integral(float_type(0),
pi<float_type>(),
tol,
[&x, &n](const float_type& t) -> float_type
{
return cos(x * sin(t) - (t * n));
});
return cos(x * sin(t) - (t * n));
}
);

return integration_result / pi<float_type>();
return integration_result / my_pi;
}
} // namespace example011_trig

Expand All @@ -392,9 +408,9 @@ auto ::math::wide_decimal::example011_trig_trapezoid_integral() -> bool
{
using example011_trig::dec51_t;

const dec51_t j2 = example011_trig::cyl_bessel_j(2U, dec51_t(123U) / 100U);
const dec51_t j3 = example011_trig::cyl_bessel_j(3U, dec51_t(456U) / 100U);
const dec51_t j4 = example011_trig::cyl_bessel_j(4U, dec51_t(789U) / 100U);
const auto j2 = example011_trig::cyl_bessel_j(static_cast<unsigned>(UINT8_C(2)), dec51_t(static_cast<unsigned>(UINT16_C(123))) / static_cast<unsigned>(UINT8_C(100)));
const auto j3 = example011_trig::cyl_bessel_j(static_cast<unsigned>(UINT8_C(3)), dec51_t(static_cast<unsigned>(UINT16_C(456))) / static_cast<unsigned>(UINT8_C(100)));
const auto j4 = example011_trig::cyl_bessel_j(static_cast<unsigned>(UINT8_C(4)), dec51_t(static_cast<unsigned>(UINT16_C(789))) / static_cast<unsigned>(UINT8_C(100)));

using std::fabs;

Expand All @@ -409,20 +425,24 @@ auto ::math::wide_decimal::example011_trig_trapezoid_integral() -> bool

using std::fabs;

const dec51_t closeness2 = fabs(1 - (j2 / control2));
const dec51_t closeness3 = fabs(1 - (j3 / control3));
const dec51_t closeness4 = fabs(1 - (j4 / control4));
const dec51_t closeness2 = fabs(static_cast<int>(INT8_C(1)) - (j2 / control2));
const dec51_t closeness3 = fabs(static_cast<int>(INT8_C(1)) - (j3 / control3));
const dec51_t closeness4 = fabs(static_cast<int>(INT8_C(1)) - (j4 / control4));

// Ensure that the local sin() function is covered for negative argument.
using example011_trig::sin;
const dec51_t closeness_sn = fabs(1 - (sin(dec51_t(-123) / 100U) / control_sn));
const dec51_t closeness_sn = fabs(1 - (sin(dec51_t(static_cast<int>(INT16_C(-123))) / static_cast<unsigned>(UINT8_C(100))) / control_sn));

const dec51_t tol = std::numeric_limits<dec51_t>::epsilon() * static_cast<std::uint32_t>(UINT8_C(10));

const auto result_is_ok = ( (closeness2 < tol)
&& (closeness3 < tol)
&& (closeness4 < tol)
&& (closeness_sn < tol));
const auto
result_is_ok =
(
(closeness2 < tol)
&& (closeness3 < tol)
&& (closeness4 < tol)
&& (closeness_sn < tol)
);

return result_is_ok;
}
Expand Down