how do I make a portable isnan/isinf function

Viewed 42297

I've been using isinf, isnan functions on Linux platforms which worked perfectly. But this didn't work on OS-X, so I decided to use std::isinf std::isnan which works on both Linux and OS-X.

But the Intel compiler doesn't recognize it, and I guess its a bug in the intel compiler according to http://software.intel.com/en-us/forums/showthread.php?t=64188

So now I just want to avoid the hassle and define my own isinf, isnan implementation.

Does anyone know how this could be done?

edit:

I ended up doing this in my source code for making isinf/isnan working

#include <iostream>
#include <cmath>

#ifdef __INTEL_COMPILER
#include <mathimf.h>
#endif

int isnan_local(double x) { 
#ifdef __INTEL_COMPILER
  return isnan(x);
#else
  return std::isnan(x);
#endif
}

int isinf_local(double x) { 
#ifdef __INTEL_COMPILER
  return isinf(x);
#else
  return std::isinf(x);
#endif
}


int myChk(double a){
  std::cerr<<"val is: "<<a <<"\t";
  if(isnan_local(a))
    std::cerr<<"program says isnan";
  if(isinf_local(a))
    std::cerr<<"program says isinf";
  std::cerr<<"\n";
  return 0;
}

int main(){
  double a = 0;
  myChk(a);
  myChk(log(a));
  myChk(-log(a));
  myChk(0/log(a));
  myChk(log(a)/log(a));

  return 0;
}
11 Answers

I've not tried this, but I would think

int isnan(double x) { return x != x; }
int isinf(double x) { return !isnan(x) && isnan(x - x); }

would work. It feels like there should be a better way for isinf, but that should work.

Well, ideally, you'd wait until Intel fixes the bug or provides a workaround :-)

But if you want to detect NaN and Inf from IEEE754 values, map it to an integer (32 or 64 bit depending on whether it's single or double precision) and check if the exponent bits are all 1. This indicates those two cases.

You can distinguish between NaN and Inf by checking the high order bit of the mantissa. If it's 1, that's NaN otherwise Inf.

+/-Inf is dictated by the sign bit.

For single precision (32-bit values), the sign is the high-order bit (b31), exponent is the next eight bits (plus a 23-bit mantissa). For double precision, the sign is still the high-order bit but the exponent is eleven bits (plus 52 bits for the mantissa).

Wikipedia has all the gory details.

The following code shows you how the encoding works.

#include <stdio.h>

static void decode (char *s, double x) {
    long y = *(((long*)(&x))+1);

    printf("%08x ",y);
    if ((y & 0x7ff80000L) == 0x7ff80000L) {
        printf ("NaN  (%s)\n", s);
        return;
    }
    if ((y & 0xfff10000L) == 0x7ff00000L) {
        printf ("+Inf (%s)\n", s);
        return;
    }
    if ((y & 0xfff10000L) == 0xfff00000L) {
        printf ("-Inf (%s)\n", s);
        return;
    }
    printf ("%e (%s)\n", x, s);
}

int main (int argc, char *argv[]) {
    double dvar;

    printf ("sizeof double = %d\n", sizeof(double));
    printf ("sizeof long   = %d\n", sizeof(long));

    dvar = 1.79e308; dvar = dvar * 10000;
    decode ("too big", dvar);

    dvar = -1.79e308; dvar = dvar * 10000;
    decode ("too big and negative", dvar);

    dvar = -1.0; dvar = sqrt(dvar);
    decode ("imaginary", dvar);

    dvar = -1.79e308;
    decode ("normal", dvar);

    return 0;
}

and it outputs:

sizeof double = 8
sizeof long   = 4
7ff00000 +Inf (too big)
fff00000 -Inf (too big and negative)
fff80000 NaN  (imaginary)
ffefdcf1 -1.790000e+308 (normal)

Just keep in mind that this code (but not the method) depends a great deal on the sizes of your longs which is not overly portable. But, if you have to bit-fiddle to get the information, you've already entered that territory :-)

As an aside, I've always found Harald Schmidt's IEEE754 converter very useful for floating point analysis.

Related