Как вернуть NaN и Infinity, проверить на NaN и Infinity?

double doubleNaN = std::numeric_limits<double>::quiet_NaN();
int intInfinity = std::numeric_limits<int>::infinity();

double doubleNaN = std::numeric_limits<double>::quiet_NaN();
double intInfinity = std::numeric_limits<double>::infinity();

typedef double (*dp)(double);
double d(dp f, double x, double dx)
{
 double result, d,oldresult;
 for (result=std::numeric_limits<double>::quiet_NaN(), d=1.0; d!=0.0; dx/=2)
 {
  if (...)
  {
   return std::numeric_limits<double>::quiet_NaN();
  }
  if (...)
  {
   return std::numeric_limits<double>::infinity();
  }
  oldresult=result;
  result=(f(x+dx)-f(x-dx))/(2*dx);
  if (/*Как написать условие oldresult!=NaN?*/)
  {
   d=oldresult-result;
  }
  else
  {
   d=1.0;
  }
 }
 return result;
}

if (y==std::numeric_limits<double>::quiet_NaN())
{
 ...
}

//http://www.johndcook.com/IEEE_exceptions_in_cpp.html
 
bool IsNumber(double x) 
    {
        // This looks like it should always be true, 
        // but it's false if x is a NaN.
        return (x == x); 
    }
 
//Here DBL_MAX is a constant defined in float.h as the largest double that can be represented.
 
bool IsFiniteNumber(double x) 
    {
        return (x <= DBL_MAX && x >= -DBL_MAX); 
    }   
 
непрохождение первого теста говорит о том что x не число (например результат ln(-1))
непрохождение второго - о том, что x бесконечность (+ или -) в зависимости от знаков операндов в операции в которой получен результат. 
 
To get more detail about the type of a floating point number, there is a function _fpclass on Windows and a corresponding function fp_class_d on Linux. I have not been able to get the corresponding Linux code to work and so I'll stick to what I've tested and just talk about Windows from here on out.
The Windows function _fpclass returns one of the following values:
 
    _FPCLASS_SNAN   // signaling NaN
        _FPCLASS_QNAN   // quiet NaN
        _FPCLASS_NINF   // negative infinity
        _FPCLASS_NN     // negative normal
        _FPCLASS_ND     // negative denormal
        _FPCLASS_NZ     // -0
        _FPCLASS_PZ     // +0
        _FPCLASS_PD     // positive denormal
        _FPCLASS_PN     // positive normal
        _FPCLASS_PINF   // positive infinity  
 
The following code illustrates which kinds of operations result in which kinds of numbers. To port this code to Linux, the FPClass function would need to use fp_class_d and its corresponding constants. 
 
 
    #include <cfloat>
        #include <iostream>
        #include <sstream>
        #include <cmath>
 
        using namespace std;
 
        string FPClass(double x)
        {
            int i = _fpclass(x);
            string s;
            switch (i)
            {
            case _FPCLASS_SNAN: s = "Signaling NaN";                break;
            case _FPCLASS_QNAN: s = "Quiet NaN";                    break; 
            case _FPCLASS_NINF: s = "Negative infinity (-INF)";     break; 
            case _FPCLASS_NN:   s = "Negative normalized non-zero"; break;
            case _FPCLASS_ND:   s = "Negative denormalized";        break; 
            case _FPCLASS_NZ:   s = "Negative zero (-0)";           break; 
            case _FPCLASS_PZ:   s = "Positive 0 (+0)";              break; 
            case _FPCLASS_PD:   s = "Positive denormalized";        break; 
            case _FPCLASS_PN:   s = "Positive normalized non-zero"; break; 
            case _FPCLASS_PINF: s = "Positive infinity (+INF)";     break;
            }
            return s;
        }
 
        string HexDump(double x)
        {
            unsigned long* pu;
            pu = (unsigned long*)&x;
            ostringstream os;
            os << hex << pu[0] << " " << pu[1];
            return os.str();
        }
 
        // ----------------------------------------------------------------------------
        int main()
        {
            double x, y, z;
 
            cout << "Testing z = 1/0\n";
            // cannot set x = 1/0 directly or would produce compile error.
            x = 1.0; y = 0; z = x/y;
            cout << "z = " << x/y << "\n";
            cout << HexDump(z) << " _fpclass(z) = " << FPClass(z) << "\n";
 
            cout << "\nTesting z = -1/0\n";
            x = -1.0; y = 0; z = x/y;
            cout << "z = " << x/y << "\n";
            cout << HexDump(z) << " _fpclass(z) = " << FPClass(z) << "\n";
 
            cout << "\nTesting z = sqrt(-1)\n";
            x = -1.0;
            z = sqrt(x);
            cout << "z = " << z << "\n";
            cout << HexDump(z) << " _fpclass(z) = " << FPClass(z) << "\n";
 
            cout << "\nTesting z = log(-1)\n";
            x = -1.0;
            z = sqrt(x);
            cout << "z = " << z << "\n";
            cout << HexDump(z) << " _fpclass(z) = " << FPClass(z) << "\n";
 
            cout << "\nTesting overflow\n";
            z = DBL_MAX;
            cout << "z = DBL_MAX = " << z; 
            z *= 2.0;
            cout << "; 2z = " << z << "\n";
            cout << HexDump(z) << " _fpclass(z) = " << FPClass(z) << "\n";
 
            cout << "\nTesting denormalized underflow\n";
            z = DBL_MIN;
            cout << "z = DBL_MIN = " << z << "\n";
            cout << HexDump(z) << " _fpclass(z) = " << FPClass(z) << "\n";
            z /= pow(2.0, 52);
            cout << "z = DBL_MIN / 2^52= " << z << "\n";
            cout << HexDump(z) << " _fpclass(z) = " << FPClass(z) << "\n";
            z /= 2;
            cout << "z = DBL_MIN / 2^53= " << z << "\n";
            cout << HexDump(z) << " _fpclass(z) = " << FPClass(z) << "\n";
 
            cout << "\nTesting z = +infinity + -infinty\n";
            x = 1.0; y = 0.0; x /= y; y = -x;
            cout << x << " + " << y << " = " << z << "\n";
            cout << HexDump(z) << " _fpclass(z) = " << FPClass(z) << "\n";
 
            cout << "\nTesting z = 0 * infinity\n";
            x = 1.0; y = 0.0; x /= y; z = 0.0*x;
            cout << "x = " << x << "; z = 0*x = " << z << "\n";
            cout << HexDump(z) << " _fpclass(z) = " << FPClass(z) << "\n";
 
            cout << "\nTesting 0/0\n";
            x = 0.0; y = 0.0; z = x/y;
            cout << "z = 0/0 = " << z << "\n";
            cout << HexDump(z) << " _fpclass(z) = " << FPClass(z) << "\n";
 
            cout << "\nTesting z = infinity/infinity\n";
            x = 1.0; y = 0.0; x /= y; y = x; z = x/y;
            cout << "x = " << x << "; z = x/x = " << z << "\n";
            cout << HexDump(z) << " _fpclass(z) = " << FPClass(z) << "\n";
 
            cout << "\nTesting x fmod 0\n";
            x = 1.0; y = 0.0; z = fmod(x, y);
            cout << "fmod(" << x << ", " << y << ") = " << z << "\n";
            cout << HexDump(z) << " _fpclass(z) = " << FPClass(z) << "\n";
 
            cout << "\nTesting infinity fmod x\n";
            y = 1.0; x = 0.0; y /= x; z = fmod(y, x);
            cout << "fmod(" << y << ", " << x << ") = " << z << "\n";
            cout << HexDump(z) << " _fpclass(z) = " << FPClass(z) << "\n";
 
            cout << "\nGetting cout to print QNAN\n";
            unsigned long nan[2]={0xffffffff, 0x7fffffff};
            z = *( double* )nan;
            cout << "z = " << z << "\n";
            cout << HexDump(z) << " _fpclass(z) = " << FPClass(z) << "\n";
 
            return 0;
        }      
 
 
немного о представлении:
// numeric_limits_infinity.cpp
// compile with: /EHsc
#include <iostream>
#include <limits>
 
using namespace std;
 
int main( )
{
   cout << numeric_limits<float>::has_infinity <<endl;
   cout << numeric_limits<double>::has_infinity<<endl;
   cout << numeric_limits<long double>::has_infinity <<endl;
   cout << numeric_limits<int>::has_infinity <<endl;
   cout << numeric_limits<__int64>::has_infinity <<endl;
 
   cout << "The representation of infinity for type float is: "
        << numeric_limits<float>::infinity( ) <<endl;
   cout << "The representation of infinity for type double is: "
        << numeric_limits<double>::infinity( ) <<endl;
   cout << "The representation of infinity for type long double is: "
        << numeric_limits<long double>::infinity( ) <<endl;
}
 
и 
Константа NaN смысла не имеет: все равно сравнивать значения с ней нельзя.
 
Константы +-Infinity в C++ получить легко:
 
std::numeric_limits<double>::infinity();
-std::numeric_limits<double>::infinity();
.

-std::numeric_limits<double>::infinity();