isnan in fortran

[c.-2]

I am regularly building Octave on OSX withg gcc/gfortran 4.2.1

since this changeset:

-----------------------------------
changeset:   14392:d17237256856
user:        Michael Goffioul <[hidden email]>
date:        Thu Feb 23 09:12:47 2012 +0000
summary:     Make SLATEC-FN atanh/erfc functions more tolerant about edge cases like
-----------------------------------

I need to locally apply the following patch (which I keep in my hq patch queue), as gfortran 4.2 does not have isnan:

-----------------------------------
# HG changeset patch
# User Carlo de Falco <[hidden email]>
# Date 1334054116 -7200
# Node ID b0734d15d631064e8963ebb4fcb3be1ead97d635
# Parent  b2bf5896ab51fdd40b899d557b971fcd0ee21994
[mq]: isnan

diff -r b2bf5896ab51 -r b0734d15d631 libcruft/slatec-fn/derfc.f
--- a/libcruft/slatec-fn/derfc.f Mon Apr 09 21:11:18 2012 +0200
+++ b/libcruft/slatec-fn/derfc.f Tue Apr 10 12:35:16 2012 +0200
@@ -191,7 +191,7 @@
       ENDIF
       FIRST = .FALSE.
 C
-      IF (ISNAN(X)) THEN
+      IF (.NOT.(X.EQ.X)) THEN
          DERFC = X
          RETURN
       ENDIF
diff -r b2bf5896ab51 -r b0734d15d631 libcruft/slatec-fn/erfc.f
--- a/libcruft/slatec-fn/erfc.f Mon Apr 09 21:11:18 2012 +0200
+++ b/libcruft/slatec-fn/erfc.f Tue Apr 10 12:35:16 2012 +0200
@@ -121,7 +121,7 @@
       ENDIF
       FIRST = .FALSE.
 C
-      IF (ISNAN(X)) THEN
+      IF (.NOT.(X.EQ.X)) THEN
          ERFC = X
          RETURN
       ENDIF
-----------------------------------

Does anyone know a cleaner and portable way how these functions could be made to work with fortran compilers that do not have ISNAN?


c.

[Mike Miller-2]

On Tue, Apr 10, 2012 at 12:44:07PM +0200, c. wrote:
> I am regularly building Octave on OSX withg gcc/gfortran 4.2.1
> [...]
> I need to locally apply the following patch (which I keep in my hq
> patch queue), as gfortran 4.2 does not have isnan:
> [...]
> Does anyone know a cleaner and portable way how these functions could
> be made to work with fortran compilers that do not have ISNAN?

I run into the same problem, RHEL 5.x with gcc 4.1.2.

Your patch is simpler and more correct than mine, which was to just
revert that cset or just comment out the entire IF..ENDIF.  Unless I
happen to have the gcc44-gfortran rpm installed, which Red Hat does
provide but is not in a standard install.

+1 to replacing this with a portable F77 solution.  If there is one...

A quick search suggests (X.NE.X) may be the best alternative, at least
for IEEE 754 systems.

--
mike

[Jordi Gutiérrez Hermoso-2]

On 10 April 2012 08:28, Mike Miller <[hidden email]> wrote:

> On Tue, Apr 10, 2012 at 12:44:07PM +0200, c. wrote:
>> I am regularly building Octave on OSX withg gcc/gfortran 4.2.1
>> [...]
>> I need to locally apply the following patch (which I keep in my hq
>> patch queue), as gfortran 4.2 does not have isnan:
>> [...]
>> Does anyone know a cleaner and portable way how these functions
>> could be made to work with fortran compilers that do not have
>> ISNAN?
>
> I run into the same problem, RHEL 5.x with gcc 4.1.2.

[snip]
> A quick search suggests (X.NE.X) may be the best alternative, at
> least for IEEE 754 systems.

Can Fortran be preprocessed? If so, can we write an autoconf macro to
check for isnan and use X.NE.X if not?

- Jordi G. H.

[Michael Goffioul]

On Tue, Apr 10, 2012 at 2:03 PM, Jordi Gutiérrez Hermoso
<[hidden email]> wrote:

> On 10 April 2012 08:28, Mike Miller <[hidden email]> wrote:
>> On Tue, Apr 10, 2012 at 12:44:07PM +0200, c. wrote:
>>> I am regularly building Octave on OSX withg gcc/gfortran 4.2.1
>>> [...]
>>> I need to locally apply the following patch (which I keep in my hq
>>> patch queue), as gfortran 4.2 does not have isnan:
>>> [...]
>>> Does anyone know a cleaner and portable way how these functions
>>> could be made to work with fortran compilers that do not have
>>> ISNAN?
>>
>> I run into the same problem, RHEL 5.x with gcc 4.1.2.
> [snip]
>> A quick search suggests (X.NE.X) may be the best alternative, at
>> least for IEEE 754 systems.
>
> Can Fortran be preprocessed? If so, can we write an autoconf macro to
> check for isnan and use X.NE.X if not?

http://gcc.gnu.org/onlinedocs/gfortran/Preprocessing-Options.html

Apparently, it can be preprocessed, but it's not enabled by default
for .f files (assuming .f extension is treated differently than .F,
which is not clear from the documentation).

Michael.

[Svante Signell]

On Tue, 2012-04-10 at 14:13 +0100, Michael Goffioul wrote:
> On Tue, Apr 10, 2012 at 2:03 PM, Jordi Gutiérrez Hermoso
...
> > Can Fortran be preprocessed? If so, can we write an autoconf macro to
> > check for isnan and use X.NE.X if not?
>
> http://gcc.gnu.org/onlinedocs/gfortran/Preprocessing-Options.html
>
> Apparently, it can be preprocessed, but it's not enabled by default
> for .f files (assuming .f extension is treated differently than .F,
> which is not clear from the documentation).

From man gcc:

       file.f
       file.for
       file.ftn
Fixed form Fortran source code which should not be preprocessed.

       file.F
       file.FOR
       file.fpp
       file.FPP
       file.FTN
Fixed form Fortran source code which must be preprocessed (with
the traditional preprocessor).

[Mike Miller-2]

On Tue, Apr 10, 2012 at 9:19 AM, Svante Signell wrote:

> On Tue, 2012-04-10 at 14:13 +0100, Michael Goffioul wrote:
>> On Tue, Apr 10, 2012 at 2:03 PM, Jordi Gutiérrez Hermoso
> ...
>> > Can Fortran be preprocessed? If so, can we write an autoconf macro to
>> > check for isnan and use X.NE.X if not?
>>
>> http://gcc.gnu.org/onlinedocs/gfortran/Preprocessing-Options.html
>>
>> Apparently, it can be preprocessed, but it's not enabled by default
>> for .f files (assuming .f extension is treated differently than .F,
>> which is not clear from the documentation).
>
> From man gcc:
>
>       file.f
>       file.for
>       file.ftn
> Fixed form Fortran source code which should not be preprocessed.
>
>       file.F
>       file.FOR
>       file.fpp
>       file.FPP
>       file.FTN
> Fixed form Fortran source code which must be preprocessed (with
> the traditional preprocessor).

Yes, and looks like automake will pass the right options to preprocess
Fortran, but only if the extension is .F, the others listed above are
possibly ignored by automake, even if gcc will handle them.

http://www.gnu.org/software/automake/manual/automake.html#Preprocessing-Fortran-77

Should work, I'll come up with something if no one beats me to it.

--
mike

[John W. Eaton]

On 10-Apr-2012, Mike Miller wrote:

| On Tue, Apr 10, 2012 at 9:19 AM, Svante Signell wrote:
| > On Tue, 2012-04-10 at 14:13 +0100, Michael Goffioul wrote:
| >> On Tue, Apr 10, 2012 at 2:03 PM, Jordi Gutiérrez Hermoso
| > ...
| >> > Can Fortran be preprocessed? If so, can we write an autoconf macro to
| >> > check for isnan and use X.NE.X if not?
| >>
| >> http://gcc.gnu.org/onlinedocs/gfortran/Preprocessing-Options.html
| >>
| >> Apparently, it can be preprocessed, but it's not enabled by default
| >> for .f files (assuming .f extension is treated differently than .F,
| >> which is not clear from the documentation).
| >
| > From man gcc:
| >
| >       file.f
| >       file.for
| >       file.ftn
| > Fixed form Fortran source code which should not be preprocessed.
| >
| >       file.F
| >       file.FOR
| >       file.fpp
| >       file.FPP
| >       file.FTN
| > Fixed form Fortran source code which must be preprocessed (with
| > the traditional preprocessor).
|
| Yes, and looks like automake will pass the right options to preprocess
| Fortran, but only if the extension is .F, the others listed above are
| possibly ignored by automake, even if gcc will handle them.
|
| http://www.gnu.org/software/automake/manual/automake.html#Preprocessing-Fortran-77
|
| Should work, I'll come up with something if no one beats me to it.

Gfortran is not the only compiler we might use to compile the Fortran
bits of Octave.

I think it would be better to rename the files that need this
treatment to FILE.in.f and then use sed to do the transformation to
generate FILE.f in the build tree.

jwe

[Mike Miller-2]

On Tue, Apr 10, 2012 at 11:19 AM, John W. Eaton <[hidden email]> wrote:

> On 10-Apr-2012, Mike Miller wrote:
>
> | On Tue, Apr 10, 2012 at 9:19 AM, Svante Signell wrote:
> | > On Tue, 2012-04-10 at 14:13 +0100, Michael Goffioul wrote:
> | >> On Tue, Apr 10, 2012 at 2:03 PM, Jordi Gutiérrez Hermoso
> | > ...
> | >> > Can Fortran be preprocessed? If so, can we write an autoconf macro to
> | >> > check for isnan and use X.NE.X if not?
> | >>
> | >> http://gcc.gnu.org/onlinedocs/gfortran/Preprocessing-Options.html
> | >>
> | >> Apparently, it can be preprocessed, but it's not enabled by default
> | >> for .f files (assuming .f extension is treated differently than .F,
> | >> which is not clear from the documentation).
> | >
> | > From man gcc:
> | >
> | >       file.f
> | >       file.for
> | >       file.ftn
> | > Fixed form Fortran source code which should not be preprocessed.
> | >
> | >       file.F
> | >       file.FOR
> | >       file.fpp
> | >       file.FPP
> | >       file.FTN
> | > Fixed form Fortran source code which must be preprocessed (with
> | > the traditional preprocessor).
> |
> | Yes, and looks like automake will pass the right options to preprocess
> | Fortran, but only if the extension is .F, the others listed above are
> | possibly ignored by automake, even if gcc will handle them.
> |
> | http://www.gnu.org/software/automake/manual/automake.html#Preprocessing-Fortran-77
> |
> | Should work, I'll come up with something if no one beats me to it.
>
> Gfortran is not the only compiler we might use to compile the Fortran
> bits of Octave.

Completely agree.  Preprocessing is by no means standard F77.
Automake has some comments that their support for it assumes the
compiler can do the work.  I have personally worked with g77,
gfortran, PG pgf77, and Intel ifort, and all support at least some of
the same preprocessing directives.  But, we can do the preprocessing
ahead anyway...

> I think it would be better to rename the files that need this
> treatment to FILE.in.f and then use sed to do the transformation to
> generate FILE.f in the build tree.

How about cpp instead of sed, so we can still have a
#if..#else..#endif syntax?  There will need to be a custom list of -D
options or a custom include file, cannot just use the same config.h,
but it's clearer and more scalable if other such feature tests are
needed in the future.

--
mike

[John W. Eaton]

On 10-Apr-2012, Mike Miller wrote:

| How about cpp instead of sed, so we can still have a
| #if..#else..#endif syntax?  There will need to be a custom list of -D
| options or a custom include file, cannot just use the same config.h,
| but it's clearer and more scalable if other such feature tests are
| needed in the future.

I'm not sure I like the idea of processing non-C syntax with cpp.

I'd rather use sed or awk for this purpose.

jwe

[c.-2]

On 10 Apr 2012, at 19:00, [hidden email] wrote:

> I think it would be better to rename the files that need this
> treatment to FILE.in.f and then use sed to do the transformation to
> generate FILE.f in the build tree.
>
> jwe

I am preparing a patch to do that.
c.

[John W. Eaton]

On 11-Apr-2012, c. wrote:

|
| On 10 Apr 2012, at 19:00, [hidden email] wrote:
|
| > I think it would be better to rename the files that need this
| > treatment to FILE.in.f and then use sed to do the transformation to
| > generate FILE.f in the build tree.
| >
| > jwe
|
| I am preparing a patch to do that.

After giving it some more thought, it might be better in this specific
case to have a configure test that checks whether isnan is available
in Fortran and if not, adds misc/isnan.f to the list of functions to
compile.  The misc/isnan.f function could be

      logical function isnan (x)
      double precision x
      if (x != x)
        isnan = .false.
      else
        isnan = .true.
      endif
      return
      end

jwe

[John W. Eaton]

On 11-Apr-2012, John W. Eaton wrote:

| On 11-Apr-2012, c. wrote:
|
| |
| | On 10 Apr 2012, at 19:00, [hidden email] wrote:
| |
| | > I think it would be better to rename the files that need this
| | > treatment to FILE.in.f and then use sed to do the transformation to
| | > generate FILE.f in the build tree.
| | >
| | > jwe
| |
| | I am preparing a patch to do that.
|
| After giving it some more thought, it might be better in this specific
| case to have a configure test that checks whether isnan is available
| in Fortran and if not, adds misc/isnan.f to the list of functions to
| compile.  The misc/isnan.f function could be
|
|       logical function isnan (x)
|       double precision x
|       if (x != x)
|         isnan = .false.
|       else
|         isnan = .true.
|       endif
|       return
|       end

Uh, except that I have the test and some syntax wrong.  How about

       logical function isnan (x)
       double precision x
       isnan = x .ne. x
       return
       end

jwe

[Michael Goffioul]

On Wed, Apr 11, 2012 at 12:24 PM, John W. Eaton <[hidden email]> wrote:

> On 11-Apr-2012, John W. Eaton wrote:
>
> | On 11-Apr-2012, c. wrote:
> |
> | |
> | | On 10 Apr 2012, at 19:00, [hidden email] wrote:
> | |
> | | > I think it would be better to rename the files that need this
> | | > treatment to FILE.in.f and then use sed to do the transformation to
> | | > generate FILE.f in the build tree.
> | | >
> | | > jwe
> | |
> | | I am preparing a patch to do that.
> |
> | After giving it some more thought, it might be better in this specific
> | case to have a configure test that checks whether isnan is available
> | in Fortran and if not, adds misc/isnan.f to the list of functions to
> | compile.  The misc/isnan.f function could be
> |
> |       logical function isnan (x)
> |       double precision x
> |       if (x != x)
> |         isnan = .false.
> |       else
> |         isnan = .true.
> |       endif
> |       return
> |       end
>
> Uh, except that I have the test and some syntax wrong.  How about
>
>       logical function isnan (x)
>       double precision x
>       isnan = x .ne. x
>       return
>       end

Sounds like a good plan. Would we also need a single precision version?

Michael.

[c.-2]

On 11 Apr 2012, at 13:24, John W. Eaton wrote:

> Uh, except that I have the test and some syntax wrong.  How about
>
>       logical function isnan (x)
>       double precision x
>       isnan = x .ne. x
>       return
>       end
>
> jwe


OK, so what about the configure test below?
How do I get to add misc/isnan.f to the list of files to compile if the test fails?
c.


diff -r b2bf5896ab51 configure.ac
--- a/configure.ac Mon Apr 09 21:11:18 2012 +0200
+++ b/configure.ac Wed Apr 11 13:32:32 2012 +0200
@@ -593,6 +593,8 @@
   AC_MSG_ERROR([in order to build octave, you must have a compatible Fortran compiler or wrapper script for f2c that functions as a Fortran compiler installed and in your path.  See the file INSTALL for more information.])
 fi
 
+OCTAVE_CHECK_FORTRAN_HAVE_ISNAN
+
 OCTAVE_CHECK_FORTRAN_INTEGER_SIZE
 if test "x$octave_cv_fortran_integer_size" = xno; then
   if $USE_64_BIT_IDX_T; then
diff -r b2bf5896ab51 m4/acinclude.m4
--- a/m4/acinclude.m4 Mon Apr 09 21:11:18 2012 +0200
+++ b/m4/acinclude.m4 Wed Apr 11 13:32:32 2012 +0200
@@ -329,6 +329,29 @@
   fi
 ])
 dnl
+dnl Check to see whether Fortran compiler has the intrinsic function ISNAN.
+dnl
+AC_DEFUN([OCTAVE_CHECK_FORTRAN_HAVE_ISNAN], [
+  AC_LANG_PUSH(Fortran 77)
+dnl        AC_CACHE_CHECK
+dnl          (
+dnl            [whether $F77 has the intrinsic function ISNAN]
+dnl            ,
+dnl            [octave_cv_fortran_have_isnan]
+dnl            ,
+dnl            [
+              AC_COMPILE_IFELSE(
+[      program foo
+      double precision x
+      logical y
+      y = isnan(x)
+      end], [AC_DEFINE([HAVE_FORTRAN_ISNAN], [1], [Define to 1 if your fortran compiler has ISNAN.])],
+                )
+dnl            ]
+dnl          )
+  AC_LANG_POP(Fortran 77)      
+])
+dnl
 dnl Check to see whether the default Fortran INTEGER is 64 bits wide.
 dnl
 AC_DEFUN([OCTAVE_CHECK_FORTRAN_INTEGER_SIZE], [

[c.-2]

On 11 Apr 2012, at 13:35, Michael Goffioul wrote:

> Sounds like a good plan. Would we also need a single precision version?
There is currently no fortran function using ISNAN for single precision in libcruft.
 
> Michael.
c.

[Michael Goffioul]

On Wed, Apr 11, 2012 at 12:39 PM, c. <[hidden email]> wrote:
>
> On 11 Apr 2012, at 13:35, Michael Goffioul wrote:
>
>> Sounds like a good plan. Would we also need a single precision version?
> There is currently no fortran function using ISNAN for single precision in libcruft.

What about
http://hg.savannah.gnu.org/hgweb/octave/file/b8edefd8c111/libcruft/slatec-fn/erfc.f

Isn't that single precision?

Michael.

[Mike Miller-2]

On Wed, Apr 11, 2012 at 07:24:18AM -0400, John W. Eaton wrote:

> On 11-Apr-2012, John W. Eaton wrote:
>
> | On 11-Apr-2012, c. wrote:
> |
> | |
> | | On 10 Apr 2012, at 19:00, [hidden email] wrote:
> | |
> | | > I think it would be better to rename the files that need this
> | | > treatment to FILE.in.f and then use sed to do the transformation to
> | | > generate FILE.f in the build tree.
> | | >
> | | > jwe
> | |
> | | I am preparing a patch to do that.
> |
> | After giving it some more thought, it might be better in this specific
> | case to have a configure test that checks whether isnan is available
> | in Fortran and if not, adds misc/isnan.f to the list of functions to
> | compile.  The misc/isnan.f function could be
> |
> |       logical function isnan (x)
> |       double precision x
> |       if (x != x)
> |         isnan = .false.
> |       else
> |         isnan = .true.
> |       endif
> |       return
> |       end
>
> Uh, except that I have the test and some syntax wrong.  How about
>
>        logical function isnan (x)
>        double precision x
>        isnan = x .ne. x
>        return
>        end

That would be a lot clearer.  So is a single-precision NaN guaranteed to
cast to a double-precision NaN value?  Or would we need isnan and
disnan?

      logical function isnan(x)
      isnan = x .ne. x
      return
      end

      logical function disnan(x)
      double precision x
      isnan = x .ne. x
      return
      end

Both erfc and derfc will need type declarations for "logical isnan" to
call this function.  Will that break the builtin isnan if the compiler
does provide it?
--
mike

[c.-2]

On 11 Apr 2012, at 13:42, Michael Goffioul wrote:

> What about
> http://hg.savannah.gnu.org/hgweb/octave/file/b8edefd8c111/libcruft/slatec-fn/erfc.f
>
> Isn't that single precision?

ouch, I'm not that proficient in F77 ...
actually there is no thype declaration for the input parameter X
I guess that means implicit typing is used, is X implicitely
typed to REAL or DOUBLE?

If that is REAL then I think we need ONLY single precision ISNAN ...

> Michael.
c.

[John W. Eaton]

On 11-Apr-2012, Mike Miller wrote:

| That would be a lot clearer.  So is a single-precision NaN guaranteed to
| cast to a double-precision NaN value?  Or would we need isnan and
| disnan?
|
|       logical function isnan(x)
|       isnan = x .ne. x
|       return
|       end
|
|       logical function disnan(x)
|       double precision x
|       isnan = x .ne. x
|       return
|       end
|
| Both erfc and derfc will need type declarations for "logical isnan" to
| call this function.  Will that break the builtin isnan if the compiler
| does provide it?

OK, that's a problem.  With the intrinsic ISNAN, the compiler chooses
the right version based on the type of the argument.  But with Fortran
77 and our own function there is no overloading, so we would need to
use different names.

jwe

[c.-2]

On 11 Apr 2012, at 14:09, John W. Eaton wrote:

> OK, that's a problem.  With the intrinsic ISNAN, the compiler chooses
> the right version based on the type of the argument.  But with Fortran
> 77 and our own function there is no overloading, so we would need to
> use different names.
>
> jwe

In that case I would rather go back to the sed or cpp option
as (X.NE.X) works with both single and double precision.
c.