Function in stdlib to get NaN

[jvdp1]

I would suggest we provide a function get_nan and internally implement it by any of the approaches that were discussed (ieee_value, sqrt or huge). That way we have just one simple place to modify, and the rest of stdlib uses get_nan and thus does not have to be modified. We should discuss a good name for such a function.

Originally posted by @certik in #128 (comment)

There was some discussion in #128 about NaN and how to generate it (not all compilers support ieee_arithmetic and its function ieee_value. Therefore, if could be nice to have such a function in stdlib that would work for all compilers.
In #128, the following solutions were proposed:

• ieee_value (from ieee_arithmetic.f90)

• sqrt of a negative number (as implemented in gcc ieee_arithmetic.f90 )

• huge (which does not return NaN??) @nshaffer

• Internal I/O to generate NaN (@urbanjost )

• transfer method (@leonfolks )

Possible name for this function: get_nan

Hopefully I didn't miss a proposition.

[fiolj]

I agree with the proposition, and would like to make it available such a function rather soon, even if the implementation will change later.
Having a function like "get_nan()" and/or a comparison function "is_nan(x)" I think would suffice for most use-cases.

I do not want to "hack" the discussion but I would suggest to have an unified discussion for the three quantities:
NaN, plus_infinite, minus_infinite

Since, they will be "special" values I would prefer an standard way to handling them rather than having each routine developer handling it individually.

[jvdp1]

I agree with the proposition, and would like to make it available such a function rather soon, even if the implementation will change later.

I agree.

Having a function like "get_nan()" and/or a comparison function "is_nan(x)" I think would suffice for most use-cases.

isnan is already supported by GFortran and Ifort, but I am not sure it is Standard Fortran, and it is only for IEEE values. So it would be indeed a nice addition.

I do not want to "hack" the discussion but I would suggest to have an unified discussion for the three quantities:
NaN, plus_infinite, minus_infinite

Good point. There are already some comments about infinity in #118 .

[nshaffer]

First, let me clarify that I suggested -huge(x) as a return value for a specific edge case for mean in order to avoid returning NaN. So we need not consider it here.

I want to emphasize that the standard is written in a way that does not assume the existence of NaN and positive/negative infinity except for IEEE real kinds. That means the sqrt, transfer, and write approaches all involve undefined behavior for non-IEEE real kinds.

So there is a fundamental friction between practical concerns (we want to be able to produce and test for NaN) and a basic concept of the language (real semantics are independent of the underlying floating-point model).

With that in mind, I propose that we have a pre-processor test for whether the ieee_arithmetic module exists. If it does, for each real kind, test if it is an IEEE kind. If it is, define NaN and +/- infinity using the IEEE intrinsic procedures. If ieee_arithmetic does not exist or if a specific kind is not IEEE-conforming, then use a fallback approach.

Of the fallback approaches discussed so far, I like the internal write best. It directly asks for a "NaN", and if the compiler knows what that means, you get something meaningful, even if it is not necessarily going to follow IEEE NaN semantics.

[certik]

Let's brainstorm the name of this function. The internal implementation of it can be decided upon later, and also changed easily.

Here are some possible names that come to my mind:

• get_nan()
• getnan()
• nan()
• ieee_nan()

Using the naming convention of "two syllables do not need a _", the getnan() seems the most natural. And also consistent with isnan().

So I vote for getnan().

Update: based on the discussion below, I now think nan() is a better choice.

[nshaffer]

Assuming we want to be able to produce NaNs of various kinds, I suggest the following API

pure elemental function nan_<kind>(x)
    real(<kind>), intent(in) :: x
    real(<kind>) :: nan_<kind>
end function nan_<kind>

where the different <kind> implementations are lumped under a generic name nan.

This approach is modeled on the intrinsic new_line(c), where the kind of the argument is used to determine the result kind. This is necessary because generic interfaces can't disambiguate based on the result alone. I prefer the name nan to reflect that it's a glorified constant. Analagously, we have new_line, not get_new_line.

I think there is no question about the API for isnan, but for completeness

pure elemental function isnan_<kind>(x)
    real(<kind>), intent(in) :: x
    logical :: isnan_<kind>
end function isnan_<kind>

The only challenge is having to be careful about playing nicely with the pre-existing (non-standard) isnan functions provided by many compilers. Presumably, we will deal with this with preprocessor flags and logic, etc.

[jvdp1]

Assuming we want to be able to produce NaNs of various kinds, I suggest the following API

pure elemental function nan_<kind>(x)
    real(<kind>), intent(in) :: x
    real(<kind>) :: nan_<kind>
end function nan_<kind>

where the different <kind> implementations are lumped under a generic name nan.

It sounds good to me. It is a similar API to huge, tiny, epsilon,...

I prefer the name nan to reflect that it's a glorified constant. Analagously, we have new_line, not get_new_line.

Based on @nshaffer reason, I vote for nan too (similarly, we have huge (not gethuge), tiny, epsilon, ...).

The only challenge is having to be careful about playing nicely with the pre-existing (non-standard) isnan functions provided by many compilers. Presumably, we will deal with this with preprocessor flags and logic, etc.

As fas as I know, when implemented, isnan supports only IEEE NaN. So it could be checked with CMake, as already done for error stop

[fiolj]

Yes, nan(), with the proposed interface sounds good to me too.
A very minor point on isnan() is that it probably should also work with complex number arguments, so the user does not have to make it manually every time.

[milancurcic]

I'm fine with either nan() or getnan() with slight preference for nan().

If we do choose to go with nan() on the basis of it being a glorified constant, I argue that we also add this specific implementation:

pure elemental function nan_default()
    real(real32) :: nan_default
end function nan_default

which would be also overloaded by the generic nan(), and which returns a single-precision real NaN if argument is omitted. Now you could really use it like a constant because Fortran allows you to call functions without (), i.e. just nan. Maybe you shouldn't though, to be consistent with the style elsewhere in the code.

Nevertheless, for people like me who almost exclusively work with single precision, nan() is nicer to use than nan(1._real32).

[certik]

I am also fine with nan and the above design.

[jvdp1]

It sounds like we have a consensus on the name and API. Based on this thread:

nan

Description

nan(x) is an elemental function that returns a NaN value of the same type as x. If x is omitted, the default type is real(sp).

Syntax

result = nan(x)

Argument

x (optional) : Shall be a scalar or an array of type integer, real or complex.

Return value

If IEEE is supported, the result is IEEE NaN and of the same type as x.
If IEEE is not supported, the result is NaN and of the same type as x.

To progress a bit on this function:

1. Note: for huge, espilon,..., the argument is not otpional. This is not a problem for me, but it might be a problem for the standard.?.

2. Implementation

• IEEE support: ieee_arithmetic.f90, ieee_value
• No IEEE support: internal I/O or transfer approach or other approaches?

The choice could be made at the compilation through CMake (and for isnan too).

[milancurcic]

What module does nan() belong to? I think the discussion first started around the stdlib_experimental_constants, but I'm not sure if that's the best place.

[nshaffer]

I don't see an issue with allowing nan to be called without arguments. However, it should return a default real, not real(sp). Our API should not make reference to specific kind parameters where it can be avoided.

I don't think the standard has any qualms with functions without arguments, see command_argument_count. It make sense for epsilon et al. to require an argument because those functions query the floating point model that a particular kind implements. The fact that new_line also requires an argument is a flaw in my opinion (though I'm curious what the standards lawyers on CLF have to say about that).

As for which module to put this in, I think stdlib_experimental_constants is OK, but I'd also support lumping NaN and infinities together in a new module stdlib_experimental_naninf.

[jvdp1]

I don't see an issue with allowing nan to be called without arguments. However, it should return a default real, not real(sp). Our API should not make reference to specific kind parameters where it can be avoided.

I see your point and agree with you for this function.

As for which module to put this in, I think stdlib_experimental_constants is OK, but I'd also support lumping NaN and infinities together in a new module stdlib_experimental_naninf.

Both are fine for me, with a preference for stdlib_experimental_naninf. E.g., do we want functions as isnan in stdlib_experimental_constants?

[nshaffer]

So this is an irrelevant detail, but I thought I'd share my findings: if nan were an intrinsic function it would not be allowed to have interfaces both with and without arguments. This is because it would be classified as an inquiry function, and inquiry functions have at least one argument. It would be an inquiry function because it depends only on properties of its argument, not the value.

I think we are not holding stdlib functions and subroutines to the same level of rigor as intrinsic ones, so the point is moot, but that's the answer.

[jvdp1]

Question:

Which compilers do not support IEEE NaN? (to test a few things)

http://fortranwiki.org/fortran/show/Fortran+2003+status

[Sideboard]

Is this being worked on?