On Wednesday, 18 December 2019 at 00:03:14 UTC, Ola Fosheim Grøstad wrote:
> On Tuesday, 17 December 2019 at 23:29:53 UTC, Timon Gehr wrote:
>> ...
>> on x86. For instance, if I added a single special case for std::pow(0.0,0.0) to a standards-compliant C++17 implementation for x86-64 with floating-point support, which values could I return without breaking C++17 standard compliance?)
>
> Whatever you like.  It is implementation defined.  That does not mean it is encouraged to return something random.
> ...

"If a domain error occurs, an implementation-defined value is returned (NaN where supported)"

I.e., what you are saying is that even if the implementation supports NaN, it may return non-NaN, the above statement notwithstanding?

> According to the standard x^y  is defined as:
>
> exp(y * log(x))
> ...

Well, that's pretty lazy. Also, it can't be true simultaneously with your claim that pow(0.0,0.0) can be modified to return _anything_, as it would then need to be consistent with exp(0.0*log(0.0)).

Also:

$ cat test.cpp
#include <cmath>
#include <iostream>
using namespace std;

int main(){
	cout<<pow(0.0,0.0)<<endl;      // 1
	cout<<exp(0.0*log(0.0))<<endl; // -nan
	double x=328.78732, y=36.3; // (random values I entered)
	cout<<(pow(x,y)==exp(y*log(x)))<<endl; // 0
}

$ g++ -std=c++11 -m64 -pedantic test.cpp && ./a.out
1
-nan
0

I may soon just go back to ignoring all your posts (like Walter also does).

On 18.12.19 01:37, Timon Gehr wrote:
> ...
> Also:
> 
> ...

Oh, and not to forget, of course exp((1.0/0.0)*log(1.0)) is NaN while pow(1.0,1.0/0.0) is 1, also invalidating the claim that allowing the implementation of pow(x,y) as exp(y*log(x)) was a goal somehow aided by treating pow(0.0,0.0) and pow(1.0,1.0/0.0) differently.

On 18.12.19 01:37, Timon Gehr wrote:
> On Wednesday, 18 December 2019 at 00:03:14 UTC, Ola Fosheim Grøstad wrote:
>> On Tuesday, 17 December 2019 at 23:29:53 UTC, Timon Gehr wrote:
>>> ...
>>> on x86. For instance, if I added a single special case for std::pow(0.0,0.0) to a standards-compliant C++17 implementation for x86-64 with floating-point support, which values could I return without breaking C++17 standard compliance?)
>>
>> Whatever you like.  It is implementation defined.  That does not mean it is encouraged to return something random.
>> ...
> 
> "If a domain error occurs, an implementation-defined value is returned (NaN where supported)"
> 
> I.e., what you are saying is that even if the implementation supports NaN, it may return non-NaN, the above statement notwithstanding?

Perhaps what you mean to say is that the C++ standard is understood to be so lax that it doesn't actually define the expected result of pow for anything but the listed special cases, such that pedantically speaking, pow could return NaN (or, usually, any other value) for all other pairs of arguments (usually, without raising a domain error)?

The webpage says that the function raises the first argument to the power of the second. For floating point, this usually means it returns the correct result rounded according to the current rounding mode. However, if it is indeed true that in the context of the C++ standard, this instead means absolutely nothing, this would successfully refute my claim that the webpage (means to) state(s) _precisely_ that pow(0.0,0.0) may return 1 or NaN after you claimed that the webpage does not say that 1 or NaN are both allowed values. It can't be true that the standard does not allow one of those two values, as NaN is explicitly allowed and actual implementations return 1.

In this case, pow(0.0,0.0) being unspecified would be exactly as significant as pow(2.0,2.0) being unspecified, and it would have exactly as much bearing on the topic of this thread. The Wikipedia article could then perhaps also be updated to explain that pow being unspecified is something that holds for most arguments, including (0.0,0.0).

On 18.12.19 03:14, Timon Gehr wrote:
> 
> Perhaps what you mean to say is that the C++ standard is understood to be so lax that it doesn't actually define the expected result of pow for anything but the listed special cases, such that pedantically speaking, pow could return NaN (or, usually, any other value) for all other pairs of arguments (usually, without raising a domain error)?

Reviewing this page, this does not appear to be the case either:
https://en.cppreference.com/w/cpp/numeric/fenv/FE_round

So I guess I still don't understand why you think an implementation could return an arbitrary value.

On 18.12.19 03:30, Timon Gehr wrote:
> On 18.12.19 03:14, Timon Gehr wrote:
>>
>> Perhaps what you mean to say is that the C++ standard is understood to be so lax that it doesn't actually define the expected result of pow for anything but the listed special cases, such that pedantically speaking, pow could return NaN (or, usually, any other value) for all other pairs of arguments (usually, without raising a domain error)?
> 
> Reviewing this page, this does not appear to be the case either:
> https://en.cppreference.com/w/cpp/numeric/fenv/FE_round
> 
> So I guess I still don't understand why you think an implementation could return an arbitrary value.

The following simple test that would have been able to refute my interpretation of the standard failed to do so:

#include <cmath>
#include <cfenv>
#include <iostream>
using namespace std;

int main(){
    #pragma STDC FENV_ACCESS ON
    double x=328.78732,y=36.3;
    fesetround(FE_DOWNWARD);
    double r1=pow(x,y);
    fesetround(FE_UPWARD);
    double r2=pow(x,y);
    cout<<*(unsigned long long*)&r1<<endl; // 5973659313751886762
    cout<<*(unsigned long long*)&r2<<endl; // 5973659313751886763
}

(Of course, this is not by itself enough to show that my interpretation is right.)

On 18.12.19 03:44, Timon Gehr wrote:
> On 18.12.19 03:30, Timon Gehr wrote:
>> On 18.12.19 03:14, Timon Gehr wrote:
>>>
>>> Perhaps what you mean to say is that the C++ standard is understood to be so lax that it doesn't actually define the expected result of pow for anything but the listed special cases, such that pedantically speaking, pow could return NaN (or, usually, any other value) for all other pairs of arguments (usually, without raising a domain error)?
>>
>> Reviewing this page, this does not appear to be the case either:
>> https://en.cppreference.com/w/cpp/numeric/fenv/FE_round
>>
>> So I guess I still don't understand why you think an implementation could return an arbitrary value.
> 
> The following simple test that would have been able to refute my interpretation of the standard failed to do so:
> 
> ...

The following D code shows that Phobos's floating point pow is a worse implementation than the one in glibc++:

import std.math, std.stdio;
void main(){
    FloatingPointControl fpctrl;
    double x=328.78732,y=36.2;
    fpctrl.rounding = FloatingPointControl.roundDown;
    double r1=x^^y;
    fpctrl.rounding = FloatingPointControl.roundUp;
    double r2=x^^y;
    writeln(*cast(ulong*)&r1); // 5969924476430611442
    writeln(*cast(ulong*)&r2); // 5969924476430611444
}

With glibc++, the two values differ by a single ulp.

On 18.12.19 03:58, Timon Gehr wrote:
> On 18.12.19 03:44, Timon Gehr wrote:
>> On 18.12.19 03:30, Timon Gehr wrote:
>>> On 18.12.19 03:14, Timon Gehr wrote:
>>>>
>>>> Perhaps what you mean to say is that the C++ standard is understood to be so lax that it doesn't actually define the expected result of pow for anything but the listed special cases, such that pedantically speaking, pow could return NaN (or, usually, any other value) for all other pairs of arguments (usually, without raising a domain error)?
>>>
>>> Reviewing this page, this does not appear to be the case either:
>>> https://en.cppreference.com/w/cpp/numeric/fenv/FE_round
>>>
>>> So I guess I still don't understand why you think an implementation could return an arbitrary value.
>>
>> The following simple test that would have been able to refute my interpretation of the standard failed to do so:
>>
>> ...

Ok, I found a case where glibc++ computes a wrong result:

#include <cmath>
#include <cfenv>
#include <iostream>
using namespace std;

int main(){
    #pragma STDC FENV_ACCESS ON
    double x=193513.887169782;
    double y=44414.97148164646;
    fesetround(FE_DOWNWARD);
    double r1=atan2(y,x);
    fesetround(FE_UPWARD);
    double r2=atan2(y,x);
    cout<<*(unsigned long long*)&r1<<endl; // 4597296506280443981
    cout<<*(unsigned long long*)&r2<<endl; // 4597296506280443981
}

If I use `long double` instead for intermediate calculations, the upper bound is the correct 4597296506280443982. So I guess the C++ standard (like IEEE floating point) does not require exact rounding for some transcendental functions, most likely including pow. Unfortunately, I haven't found any details about precision requirements for C++ floating point library functions using a cursory Google search, so it may indeed be the case that there are absolutely none. Do you have any source that says as much?

On 18.12.19 04:29, Timon Gehr wrote:
> So I guess the C++ standard (like IEEE floating point) does not require exact rounding for some transcendental functions, most likely including pow.

It seems IEEE 754 recommends exactly rounded `pow`, but it is not required for conformance (but `pow(0.0,0.0)==1.0` is required).

Also, I finally found this: https://en.cppreference.com/w/cpp/numeric/math/sqrt

"Notes
std::sqrt is required by the IEEE standard to be exact. The only other operations required to be exact are the arithmetic operators and the function std::fma. After rounding to the return type (using default rounding mode), the result of std::sqrt is indistinguishable from the infinitely precise result. In other words, the error is less than 0.5 ulp. Other functions, including std::pow, are not so constrained."

What a weird place to hide this information. I guess I was wrong about C++11 `pow(0.0,0.0)` being _required_ to return either 1.0 or NaN. Sorry about that. (However, any implementation that chooses to implement `pow` conforming to IEEE 754 will return `1.0`.)

Unfortunately, none of this is actually relevant for the case of D's `pow(int,int)`.

On 18.12.19 01:37, Timon Gehr wrote:
> 
>> According to the standard x^y  is defined as:
>>
>> exp(y * log(x))
>> ...
> 
> Well, that's pretty lazy. Also, it can't be true simultaneously with your claim that pow(0.0,0.0) can be modified to return _anything_, as it would then need to be consistent with exp(0.0*log(0.0)).

I guess what's going on is that exp and log in your expression as it occurs in the standard are the actual mathematical functions and `pow` is defined to approximate this exact result for arguments where it is defined, while for other arguments there is an explicit definition.

Another thing I noticed is that https://en.cppreference.com/w/cpp/numeric/math/pow says:
"except where specified above, if any argument is NaN, NaN is returned"

As pow(NaN,0.0) is not "specified above", this seems to say that pow(0.0/0.0,0.0) should be NaN. However, g++ gives me 1.

I'm not sure what's going on here. I guess either the documentation is wrong or g++ violates the C++ standard in order to satisfy recommendations in IEEE 754-2008.

On 18.12.19 05:23, Timon Gehr wrote:
> 
> As pow(NaN,0.0) is not "specified above", this seems to say that pow(0.0/0.0,0.0) should be NaN. However, g++ gives me 1.

https://en.cppreference.com/w/cpp/numeric/math/pow

"pow(base, ±0) returns 1 for any base, even when base is NaN"

Rofl. Ola somehow managed to gaslight me into thinking that wasn't there after I had already used it to draw conclusions. Therefore, the page _really_ states that the result must be either 1 or NaN (in case the implementation supports NaN, otherwise it can be anything), as I originally claimed, and the only thing I was wrong about tonight was I being wrong. Sorry for the noise.

Also, that does it. Ola, you are back in my kill file.