Skip to content

Finish P0811R3 midpoint and lerp #1048

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

Merged
CaseyCarter merged 12 commits into from
Jul 30, 2020
Merged

Finish P0811R3 midpoint and lerp #1048

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
12 commits
Select commit Hold shift + click to select a range
fd1590b
Finish P0811R3 midpoint(), lerp()
statementreply Jul 16, 2020
d554f60
Update tests for lerp()
statementreply Jul 16, 2020
e3a405a
Apply clang-format
statementreply Jul 16, 2020
e5c0662
Use template variable
statementreply Jul 22, 2020
e2f612b
Use assert in constexpr test
statementreply Jul 22, 2020
033632e
Add comments for sNaN handling code
statementreply Jul 23, 2020
98b19bb
Update feature test macro
statementreply Jul 23, 2020
fbd2c52
Update tests
statementreply Jul 23, 2020
9305264
Add tests for floating point exceptions
statementreply Jul 23, 2020
c9a2ba5
Minor fixes
statementreply Jul 23, 2020
57e9826
Address review comments
statementreply Jul 29, 2020
af71038
Fix a subtle conformance issue.
StephanTLavavej Jul 30, 2020
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
108 changes: 38 additions & 70 deletions stl/inc/cmath
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -12,6 +12,10 @@
#include <cstdlib>
#include <xtr1common>

#if _HAS_CXX20
#include <xutility>
#endif // _HAS_CXX20

#pragma pack(push, _CRT_PACKING)
#pragma warning(push, _STL_WARNING_LEVEL)
#pragma warning(disable : _STL_DISABLED_WARNINGS)
Expand Down Expand Up @@ -1238,13 +1242,12 @@ _NODISCARD auto hypot(const _Ty1 _Dx, const _Ty2 _Dy, const _Ty3 _Dz) {

#if _HAS_CXX20
// FUNCTION lerp
// TRANSITION, P0553: lerp is not yet constexpr
template <class _Ty>
_NODISCARD /* constexpr */ _Ty _Common_lerp(const _Ty _ArgA, const _Ty _ArgB, const _Ty _ArgT) noexcept {
_NODISCARD constexpr _Ty _Common_lerp(const _Ty _ArgA, const _Ty _ArgB, const _Ty _ArgT) noexcept {
// on a line intersecting {(0.0, _ArgA), (1.0, _ArgB)}, return the Y value for X == _ArgT

const int _Finite_mask = (int{isfinite(_ArgA)} << 2) | (int{isfinite(_ArgB)} << 1) | int{isfinite(_ArgT)};
if (_Finite_mask == 0b111) {
const bool _T_is_finite = _STD _Is_finite(_ArgT);
if (_T_is_finite && _STD _Is_finite(_ArgA) && _STD _Is_finite(_ArgB)) {
// 99% case, put it first; this block comes from P0811R3
if ((_ArgA <= 0 && _ArgB >= 0) || (_ArgA >= 0 && _ArgB <= 0)) {
// exact, monotonic, bounded, determinate, and (for _ArgA == _ArgB == 0) consistent:
Expand Down Expand Up @@ -1272,96 +1275,61 @@ _NODISCARD /* constexpr */ _Ty _Common_lerp(const _Ty _ArgA, const _Ty _ArgB, co
return _Candidate;
}

if (isnan(_ArgA)) {
return _ArgA;
}

if (isnan(_ArgB)) {
return _ArgB;
}

if (isnan(_ArgT)) {
return _ArgT;
}

switch (_Finite_mask) {
case 0b000:
// All values are infinities
if (_ArgT >= 1) {
return _ArgB;
}

return _ArgA;
case 0b010:
case 0b100:
case 0b110:
// _ArgT is an infinity; return infinity in the "direction" of _ArgA and _ArgB
return _ArgT * (_ArgB - _ArgA);
case 0b001:
// Here _ArgA and _ArgB are infinities
if (_ArgA == _ArgB) {
// same sign, so T doesn't matter
return _ArgA;
}

// Opposite signs, choose the "infinity direction" according to T if it makes sense.
if (_ArgT <= 0) {
if (_STD is_constant_evaluated()) {
if (_STD _Is_nan(_ArgA)) {
return _ArgA;
}

if (_ArgT >= 1) {
if (_STD _Is_nan(_ArgB)) {
return _ArgB;
}

// Interpolating between infinities of opposite signs doesn't make sense, NaN
if constexpr (sizeof(_Ty) == sizeof(float)) {
return __builtin_nanf("0");
} else {
return __builtin_nan("0");
}
case 0b011:
// _ArgA is an infinity but _ArgB is not
if (_ArgT == 1) {
return _ArgB;
if (_STD _Is_nan(_ArgT)) {
return _ArgT;
}

if (_ArgT < 1) {
// towards the infinity, return it
return _ArgA;
} else {
// raise FE_INVALID if at least one of _ArgA, _ArgB, and _ArgT is signaling NaN
if (_STD _Is_nan(_ArgA) || _STD _Is_nan(_ArgB)) {
return (_ArgA + _ArgB) + _ArgT;
}

// away from the infinity
return -_ArgA;
case 0b101:
// _ArgA is finite and _ArgB is an infinity
if (_ArgT == 0) {
return _ArgA;
if (_STD _Is_nan(_ArgT)) {
return _ArgT + _ArgT;
}
}

if (_ArgT > 0) {
// toward the infinity
return _ArgB;
if (_T_is_finite) {
// _ArgT is finite, _ArgA and/or _ArgB is infinity
if (_ArgT < 0) {
// if _ArgT < 0: return infinity in the "direction" of _ArgA if that exists, NaN otherwise
return _ArgA - _ArgB;
} else if (_ArgT <= 1) {
// if _ArgT == 0: return _ArgA (infinity) if _ArgB is finite, NaN otherwise
// if 0 < _ArgT < 1: return infinity "between" _ArgA and _ArgB if that exists, NaN otherwise
// if _ArgT == 1: return _ArgB (infinity) if _ArgA is finite, NaN otherwise
return _ArgT * _ArgB + (1 - _ArgT) * _ArgA;
} else {
// if _ArgT > 1: return infinity in the "direction" of _ArgB if that exists, NaN otherwise
return _ArgB - _ArgA;
}

return -_ArgB;
case 0b111: // impossible; handled in fast path
default:
_CSTD abort();
} else {
// _ArgT is an infinity; return infinity in the "direction" of _ArgA and _ArgB if that exists, NaN otherwise
return _ArgT * (_ArgB - _ArgA);
}
}

// As of 2019-06-17 it is unclear whether the "sufficient additional overloads" clause is intended to target lerp;
// LWG-3223 is pending.

_NODISCARD /* constexpr */ inline float lerp(const float _ArgA, const float _ArgB, const float _ArgT) noexcept {
_NODISCARD constexpr inline float lerp(const float _ArgA, const float _ArgB, const float _ArgT) noexcept {
return _Common_lerp(_ArgA, _ArgB, _ArgT);
}

_NODISCARD /* constexpr */ inline double lerp(const double _ArgA, const double _ArgB, const double _ArgT) noexcept {
_NODISCARD constexpr inline double lerp(const double _ArgA, const double _ArgB, const double _ArgT) noexcept {
return _Common_lerp(_ArgA, _ArgB, _ArgT);
}

_NODISCARD /* constexpr */ inline long double lerp(
_NODISCARD constexpr inline long double lerp(
const long double _ArgA, const long double _ArgB, const long double _ArgT) noexcept {
return _Common_lerp(_ArgA, _ArgB, _ArgT);
}
Expand Down
48 changes: 24 additions & 24 deletions stl/inc/numeric
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -536,21 +536,16 @@ _CONSTEXPR20 void iota(_FwdIt _First, _FwdIt _Last, _Ty _Val) {

#if _HAS_CXX17
// FUNCTION TEMPLATE _Abs_u
template <class _Arithmetic>
_NODISCARD constexpr auto _Abs_u(const _Arithmetic _Val) noexcept {
template <class _Integral>
_NODISCARD constexpr auto _Abs_u(const _Integral _Val) noexcept {
// computes absolute value of _Val (converting to an unsigned integer type if necessary to avoid overflow
// representing the negation of the minimum value)
if constexpr (is_floating_point_v<_Arithmetic>) {
// TRANSITION, P0553: this mishandles NaNs
if (_Val < 0) {
return -_Val;
}
static_assert(is_integral_v<_Integral>);

return _Val;
} else if constexpr (is_signed_v<_Arithmetic>) {
using _Unsigned = make_unsigned_t<_Arithmetic>;
if constexpr (is_signed_v<_Integral>) {
using _Unsigned = make_unsigned_t<_Integral>;
if (_Val < 0) {
// note static_cast to _Unsigned such that _Arithmetic == short returns unsigned short rather than int
// note static_cast to _Unsigned such that _Integral == short returns unsigned short rather than int
return static_cast<_Unsigned>(_Unsigned{0} - static_cast<_Unsigned>(_Val));
}

Expand Down Expand Up @@ -619,9 +614,24 @@ _NODISCARD constexpr common_type_t<_Mt, _Nt> lcm(const _Mt _Mx, const _Nt _Nx) n
template <class _Ty, enable_if_t<is_arithmetic_v<_Ty> && !is_same_v<remove_cv_t<_Ty>, bool>, int> = 0>
_NODISCARD constexpr _Ty midpoint(const _Ty _Val1, const _Ty _Val2) noexcept {
if constexpr (is_floating_point_v<_Ty>) {
if (_STD is_constant_evaluated()) {
if (_STD _Is_nan(_Val1)) {
return _Val1;
}

if (_STD _Is_nan(_Val2)) {
return _Val2;
}
} else {
if (_STD _Is_nan(_Val1) || _STD _Is_nan(_Val2)) {
Copy link
Member

@BillyONeal BillyONeal Jul 22, 2020

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Before I purposely put this test below the high limit tests before the NaN tests because we expect NaNs to be uncommon; you're fixing -fassociative-math and friends but it seems like that would still be fixed if you put this in the old location?

Copy link
Contributor Author

@statementreply statementreply Jul 23, 2020 •
edited
Loading

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Performing the high limit test with <= raises undeserved FE_INVALID when the input value is NaN.

Alternative approaches:

  • Test with quiet comparison islessequal
  • Test with bit_cast tricks
  • Just ignore the issue of FE_INVALID

Copy link
Member

@BillyONeal BillyONeal Jul 23, 2020

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I see, to raise floating point exceptions. Can you add a comment to places where you do seemingly needless operations to raise the right exceptions and add tests for that? Otherwise that's likely to break. (I understand that requires /fp:except which might be annoying to set up in our tests :( )

Copy link
Member

@BillyONeal BillyONeal Jul 23, 2020

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Hmmm actually it looks like we already have an /fp:except flavor tested:

PM_CL="/EHsc /MDd /D_ITERATOR_DEBUG_LEVEL=2 /std:c++latest /fp:except /Zc:preprocessor"

so you should just be able to do that.

// raise FE_INVALID if at least one of _Val1 and _Val2 is signaling NaN
return _Val1 + _Val2;
}
}

constexpr _Ty _High_limit = (numeric_limits<_Ty>::max)() / 2;
const auto _Val1_a = _Abs_u(_Val1);
const auto _Val2_a = _Abs_u(_Val2);
const auto _Val1_a = _Float_abs(_Val1);
const auto _Val2_a = _Float_abs(_Val2);
if (_Val1_a <= _High_limit && _Val2_a <= _High_limit) {
// _Val1 and _Val2 are small enough that _Val1 + _Val2 won't overflow

Expand All @@ -637,22 +647,12 @@ _NODISCARD constexpr _Ty midpoint(const _Ty _Val1, const _Ty _Val2) noexcept {
return (_Val1 + _Val2) / 2;
}

// TRANSITION, P0553: the next two branches handle NaNs but don't produce correct behavior under /fp:fast or
// -fassociative-math
if (_Val1 != _Val1) {
return _Val1;
}

if (_Val2 != _Val2) {
return _Val2;
}

// Here at least one of {_Val1, _Val2} has large magnitude.
// Therefore, if one of the values is too small to divide by 2 exactly, the small magnitude is much less than
// one ULP of the result, so we can add it directly without the potentially inexact division by 2.

// In the default rounding mode this less than one ULP difference will always be rounded away, so under
// /fp:precise or /fp:fast we could avoid these tests if we had some means of detecting it in the caller.
// /fp:fast we could avoid these tests if we had some means of detecting it in the caller.
constexpr _Ty _Low_limit = (numeric_limits<_Ty>::min)() * 2;
if (_Val1_a < _Low_limit) {
return _Val1 + _Val2 / 2;
Expand Down
Loading