changeset 21721:bcc30b45a225

define octave_numeric_limits template methods for Inf and NaN * lo-ieee.h (octave_numeric_limits::Inf, octave_numeric_limits::NaN): New methods. * Canvas.cc, data.cc, graphics.cc, graphics.in.h, oct-handle.h, oct-stream.cc, quadcc.cc, sparse-xdiv.cc, str2double.cc, variables.h, __glpk__.cc, __init_fltk__.cc, __voronoi__.cc, octave.cc, Matrix.cc, CNDArray.cc, CSparse.cc, Sparse-C.cc, dDiagMatrix.cc, dMatrix.cc, dSparse.cc, fCMatrix.cc, fCNDArray.cc, fDiagMatrix.cc, fMatrix.cc, Faddeeva.cc, eigs-base.cc, lo-mappers.cc, lo-mappers.h, lo-specfun.cc, oct-norm.cc, oct-rand.cc, oct-rand.h, oct-spparms.cc, randgamma.cc, randpoisson.cc: Use them instead of octave_Inf, octave_Float_Inf, octave_NaN, and octave_Float_NaN.
author John W. Eaton <jwe@octave.org>
date Tue, 17 May 2016 01:11:27 -0400
parents b28c26ae737c
children fb5dd9f7d697
files libgui/graphics/Canvas.cc libinterp/corefcn/data.cc libinterp/corefcn/graphics.cc libinterp/corefcn/graphics.in.h libinterp/corefcn/oct-handle.h libinterp/corefcn/oct-stream.cc libinterp/corefcn/quadcc.cc libinterp/corefcn/sparse-xdiv.cc libinterp/corefcn/str2double.cc libinterp/corefcn/variables.h libinterp/dldfcn/__glpk__.cc libinterp/dldfcn/__init_fltk__.cc libinterp/dldfcn/__voronoi__.cc libinterp/octave.cc liboctave/array/CMatrix.cc liboctave/array/CNDArray.cc liboctave/array/CSparse.cc liboctave/array/Sparse-C.cc liboctave/array/dDiagMatrix.cc liboctave/array/dMatrix.cc liboctave/array/dSparse.cc liboctave/array/fCMatrix.cc liboctave/array/fCNDArray.cc liboctave/array/fDiagMatrix.cc liboctave/array/fMatrix.cc liboctave/cruft/Faddeeva/Faddeeva.cc liboctave/numeric/eigs-base.cc liboctave/numeric/lo-mappers.cc liboctave/numeric/lo-mappers.h liboctave/numeric/lo-specfun.cc liboctave/numeric/oct-norm.cc liboctave/numeric/oct-rand.cc liboctave/numeric/oct-rand.h liboctave/numeric/oct-spparms.cc liboctave/numeric/randgamma.cc liboctave/numeric/randpoisson.cc liboctave/util/lo-ieee.h
diffstat 37 files changed, 255 insertions(+), 232 deletions(-) [+]
line wrap: on
line diff
--- a/libgui/graphics/Canvas.cc	Mon May 16 23:20:50 2016 -0400
+++ b/libgui/graphics/Canvas.cc	Tue May 17 01:11:27 2016 -0400
@@ -615,7 +615,7 @@
         Utils::properties<figure> (figObj)
           .set_currentobject (currentObj.get_handle ().as_octave_value ());
       else
-        Utils::properties<figure> (figObj).set_currentobject (octave_NaN);
+        Utils::properties<figure> (figObj).set_currentobject (octave_numeric_limits<double>::NaN ());
 
       Figure* fig = dynamic_cast<Figure*> (Backend::toolkitObject (figObj));
 
--- a/libinterp/corefcn/data.cc	Mon May 16 23:20:50 2016 -0400
+++ b/libinterp/corefcn/data.cc	Tue May 17 01:11:27 2016 -0400
@@ -5669,7 +5669,7 @@
       if (str == "fro")
         p_arg = octave_value (2);
       else if (str == "inf")
-        p_arg = octave_Inf;
+        p_arg = octave_numeric_limits<double>::Inf ();
       else
         error ("norm: unrecognized option: %s", str.c_str ());
     }
@@ -5697,7 +5697,7 @@
       break;
 
     case sfinf:
-      retval = xnorm (x_arg, octave_Inf);
+      retval = xnorm (x_arg, octave_numeric_limits<double>::Inf ());
       break;
     }
 
--- a/libinterp/corefcn/graphics.cc	Mon May 16 23:20:50 2016 -0400
+++ b/libinterp/corefcn/graphics.cc	Tue May 17 01:11:27 2016 -0400
@@ -1503,8 +1503,8 @@
 void
 array_property::get_data_limits (void)
 {
-  xmin = xminp = octave_Inf;
-  xmax = xmaxp = -octave_Inf;
+  xmin = xminp = octave_numeric_limits<double>::Inf ();
+  xmax = xmaxp = -octave_numeric_limits<double>::Inf ();
 
   if (! data.is_empty ())
     {
@@ -1722,7 +1722,7 @@
     }
   else if (type.compare ("handle"))
     {
-      double hv = (args.length () > 0 ? args(0).double_value () : octave_NaN);
+      double hv = (args.length () > 0 ? args(0).double_value () : octave_numeric_limits<double>::NaN ());
 
       graphics_handle gh (hv);
 
@@ -2579,7 +2579,7 @@
           const std::string& pname, const graphics_handle& new_parent,
           bool adopt = true)
 {
-  graphics_handle h = octave_NaN;
+  graphics_handle h = octave_numeric_limits<double>::NaN ();
 
   double hv = ov.xdouble_value ("%s: %s must be a graphics handle",
                                who.c_str (), pname.c_str ());
@@ -2612,7 +2612,7 @@
 {
   octave_value val = xget (0, "currentfigure");
 
-  return val.is_empty () ? octave_NaN : val.double_value ();
+  return val.is_empty () ? octave_numeric_limits<double>::NaN () : val.double_value ();
 }
 
 // This function is NOT equivalent to the scripting language function gca.
@@ -2621,7 +2621,7 @@
 {
   octave_value val = xget (gcf (), "currentaxes");
 
-  return val.is_empty () ? octave_NaN : val.double_value ();
+  return val.is_empty () ? octave_numeric_limits<double>::NaN () : val.double_value ();
 }
 
 static void
@@ -2997,7 +2997,7 @@
 {
   double hp = val.xdouble_value ("set: parent must be a graphics handle");
 
-  graphics_handle new_parent = octave_NaN;
+  graphics_handle new_parent = octave_numeric_limits<double>::NaN ();
 
   if (hp == __myhandle__)
     error ("set: can not set object parent to be object itself");
@@ -6055,10 +6055,10 @@
 {
   if (tight)
     {
-      double minval = octave_Inf;
-      double maxval = -octave_Inf;
-      double min_pos = octave_Inf;
-      double max_neg = -octave_Inf;
+      double minval = octave_numeric_limits<double>::Inf ();
+      double maxval = -octave_numeric_limits<double>::Inf ();
+      double min_pos = octave_numeric_limits<double>::Inf ();
+      double max_neg = -octave_numeric_limits<double>::Inf ();
       get_children_limits (minval, maxval, min_pos, max_neg, kids, limit_type);
       if (xfinite (minval) && xfinite (maxval))
         {
@@ -6109,7 +6109,7 @@
     }
   else
     {
-      double s = -octave_Inf;
+      double s = -octave_numeric_limits<double>::Inf ();
       bool modified_limits = false;
       Matrix kids;
 
@@ -6230,8 +6230,8 @@
   graphics_xform xform = get_transform ();
 
   Matrix ext (1, 4, 0.0);
-  ext(0) = ext(1) = octave_Inf;
-  ext(2) = ext(3) = -octave_Inf;
+  ext(0) = ext(1) = octave_numeric_limits<double>::Inf ();
+  ext(2) = ext(3) = -octave_numeric_limits<double>::Inf ();
   for (int i = 0; i <= 1; i++)
     for (int j = 0; j <= 1; j++)
       for (int k = 0; k <= 1; k++)
@@ -7114,10 +7114,10 @@
 
   Matrix kids = Matrix (1, 1, h.value ());
 
-  double min_val = octave_Inf;
-  double max_val = -octave_Inf;
-  double min_pos = octave_Inf;
-  double max_neg = -octave_Inf;
+  double min_val = octave_numeric_limits<double>::Inf ();
+  double max_val = -octave_numeric_limits<double>::Inf ();
+  double min_pos = octave_numeric_limits<double>::Inf ();
+  double max_neg = -octave_numeric_limits<double>::Inf ();
 
   char update_type = 0;
 
@@ -7318,10 +7318,10 @@
 
   Matrix kids = xproperties.get_children ();
 
-  double min_val = octave_Inf;
-  double max_val = -octave_Inf;
-  double min_pos = octave_Inf;
-  double max_neg = -octave_Inf;
+  double min_val = octave_numeric_limits<double>::Inf ();
+  double max_val = -octave_numeric_limits<double>::Inf ();
+  double min_pos = octave_numeric_limits<double>::Inf ();
+  double max_neg = -octave_numeric_limits<double>::Inf ();
 
   char update_type = 0;
 
@@ -7542,16 +7542,16 @@
 
   // Get children axes limits
   Matrix kids = get_children ();
-  double minx = octave_Inf;
-  double maxx = -octave_Inf;
-  double min_pos_x = octave_Inf;
-  double max_neg_x = -octave_Inf;
+  double minx = octave_numeric_limits<double>::Inf ();
+  double maxx = -octave_numeric_limits<double>::Inf ();
+  double min_pos_x = octave_numeric_limits<double>::Inf ();
+  double max_neg_x = -octave_numeric_limits<double>::Inf ();
   get_children_limits (minx, maxx, min_pos_x, max_neg_x, kids, 'x');
 
-  double miny = octave_Inf;
-  double maxy = -octave_Inf;
-  double min_pos_y = octave_Inf;
-  double max_neg_y = -octave_Inf;
+  double miny = octave_numeric_limits<double>::Inf ();
+  double maxy = -octave_numeric_limits<double>::Inf ();
+  double min_pos_y = octave_numeric_limits<double>::Inf ();
+  double max_neg_y = -octave_numeric_limits<double>::Inf ();
   get_children_limits (miny, maxy, min_pos_y, max_neg_y, kids, 'y');
 
   xlims = do_zoom (x, factor, xlims, xscale_is ("log"));
@@ -7689,16 +7689,16 @@
 
   // Get children axes limits
   Matrix kids = get_children ();
-  double minx = octave_Inf;
-  double maxx = -octave_Inf;
-  double min_pos_x = octave_Inf;
-  double max_neg_x = -octave_Inf;
+  double minx = octave_numeric_limits<double>::Inf ();
+  double maxx = -octave_numeric_limits<double>::Inf ();
+  double min_pos_x = octave_numeric_limits<double>::Inf ();
+  double max_neg_x = -octave_numeric_limits<double>::Inf ();
   get_children_limits (minx, maxx, min_pos_x, max_neg_x, kids, 'x');
 
-  double miny = octave_Inf;
-  double maxy = -octave_Inf;
-  double min_pos_y = octave_Inf;
-  double max_neg_y = -octave_Inf;
+  double miny = octave_numeric_limits<double>::Inf ();
+  double maxy = -octave_numeric_limits<double>::Inf ();
+  double min_pos_y = octave_numeric_limits<double>::Inf ();
+  double max_neg_y = -octave_numeric_limits<double>::Inf ();
   get_children_limits (miny, maxy, min_pos_y, max_neg_y, kids, 'y');
 
   xlims = do_translate (x0, x1, xlims, xscale_is ("log"));
@@ -8435,10 +8435,10 @@
 
   Matrix kids = Matrix (1, 1, h.value ());
 
-  double min_val = octave_Inf;
-  double max_val = -octave_Inf;
-  double min_pos = octave_Inf;
-  double max_neg = -octave_Inf;
+  double min_val = octave_numeric_limits<double>::Inf ();
+  double max_val = -octave_numeric_limits<double>::Inf ();
+  double min_pos = octave_numeric_limits<double>::Inf ();
+  double max_neg = -octave_numeric_limits<double>::Inf ();
 
   Matrix limits;
   double val;
@@ -8548,10 +8548,10 @@
 
   Matrix kids = xproperties.get_children ();
 
-  double min_val = octave_Inf;
-  double max_val = -octave_Inf;
-  double min_pos = octave_Inf;
-  double max_neg = -octave_Inf;
+  double min_val = octave_numeric_limits<double>::Inf ();
+  double max_val = -octave_numeric_limits<double>::Inf ();
+  double min_pos = octave_numeric_limits<double>::Inf ();
+  double max_neg = -octave_numeric_limits<double>::Inf ();
 
   char update_type = 0;
 
@@ -10088,7 +10088,7 @@
 {
   octave_value retval;
 
-  double val = octave_NaN;
+  double val = octave_numeric_limits<double>::NaN ();
 
   octave_value_list xargs = args.splice (0, 1);
 
@@ -10173,7 +10173,7 @@
 
       octave_value_list xargs = args.splice (0, 1);
 
-      graphics_handle h = octave_NaN;
+      graphics_handle h = octave_numeric_limits<double>::NaN ();
 
       if (xisnan (val))
         {
@@ -10427,7 +10427,7 @@
   if (args.length () != 1)
     print_usage ();
 
-  graphics_handle h = octave_NaN;
+  graphics_handle h = octave_numeric_limits<double>::NaN ();
 
   const NDArray vals = args(0).xarray_value ("delete: invalid graphics object");
 
@@ -10463,7 +10463,7 @@
   if (nargin == 2)
     mode = args(1).string_value ();
 
-  graphics_handle h = octave_NaN;
+  graphics_handle h = octave_numeric_limits<double>::NaN ();
 
   double val = args(0).xdouble_value ("__go_axes_init__: invalid graphics object");
 
--- a/libinterp/corefcn/graphics.in.h	Mon May 16 23:20:50 2016 -0400
+++ b/libinterp/corefcn/graphics.in.h	Tue May 17 01:11:27 2016 -0400
@@ -1537,7 +1537,7 @@
     return *this;
   }
 
-  void invalidate (void) { current_val = octave_NaN; }
+  void invalidate (void) { current_val = octave_numeric_limits<double>::NaN (); }
 
   base_property* clone (void) const { return new handle_property (*this); }
 
@@ -4530,17 +4530,17 @@
       Matrix lim;
 
       lim = Matrix (1, 3, pos(0));
-      lim(2) = (lim(2) <= 0 ? octave_Inf : lim(2));
+      lim(2) = (lim(2) <= 0 ? octave_numeric_limits<double>::Inf () : lim(2));
       set_xlim (lim);
 
       lim = Matrix (1, 3, pos(1));
-      lim(2) = (lim(2) <= 0 ? octave_Inf : lim(2));
+      lim(2) = (lim(2) <= 0 ? octave_numeric_limits<double>::Inf () : lim(2));
       set_ylim (lim);
 
       if (pos.numel () == 3)
         {
           lim = Matrix (1, 3, pos(2));
-          lim(2) = (lim(2) <= 0 ? octave_Inf : lim(2));
+          lim(2) = (lim(2) <= 0 ? octave_numeric_limits<double>::Inf () : lim(2));
           set_zliminclude ("on");
           set_zlim (lim);
         }
--- a/libinterp/corefcn/oct-handle.h	Mon May 16 23:20:50 2016 -0400
+++ b/libinterp/corefcn/oct-handle.h	Tue May 17 01:11:27 2016 -0400
@@ -35,10 +35,10 @@
 class octave_handle
 {
 public:
-  octave_handle (void) : val (octave_NaN) { }
+  octave_handle (void) : val (octave_numeric_limits<double>::NaN ()) { }
 
   octave_handle (const octave_value& a)
-    : val (octave_NaN)
+    : val (octave_numeric_limits<double>::NaN ())
   {
     if (a.is_empty ())
       ; // do nothing
--- a/libinterp/corefcn/oct-stream.cc	Mon May 16 23:20:50 2016 -0400
+++ b/libinterp/corefcn/oct-stream.cc	Tue May 17 01:11:27 2016 -0400
@@ -2515,7 +2515,7 @@
   : who (who_arg), buf (), whitespace_table (), delim_table (),
     delims (), comment_style (), comment_len (0), comment_char (-2),
     buffer_size (0), date_locale (), inf_nan (init_inf_nan ()),
-    empty_value (octave_NaN), exp_chars ("edED"),
+    empty_value (octave_numeric_limits<double>::NaN ()), exp_chars ("edED"),
     header_lines (0), treat_as_empty (), treat_as_empty_len (0),
     whitespace (" \b\t"), eol1 ('\r'), eol2 ('\n'),
     return_on_error (1), collect_output (false),
@@ -2899,12 +2899,12 @@
       int i = lookahead (is, inf_nan, 3, false);   // false -> case insensitive
       if (i == 0)
         {
-          retval = octave_Inf;
+          retval = octave_numeric_limits<double>::Inf ();
           valid = true;
         }
       else if (i == 1)
         {
-          retval = octave_NaN;
+          retval = octave_numeric_limits<double>::NaN ();
           valid = true;
         }
     }
@@ -2915,12 +2915,12 @@
       int i = lookahead (is, inf_nan, 3, false);   // false -> case insensitive
       if (i == 0)
         {
-          retval = octave_Inf;
+          retval = octave_numeric_limits<double>::Inf ();
           valid = true;
         }
       else if (i == 1)
         {
-          retval = octave_NaN;
+          retval = octave_numeric_limits<double>::NaN ();
           valid = true;
         }
     }
@@ -2966,7 +2966,7 @@
               if (ch2 == 'f')
                 {
                   inf = true;
-                  re = (ch == '+') ? octave_Inf : -octave_Inf;
+                  re = (ch == '+') ? octave_numeric_limits<double>::Inf () : -octave_numeric_limits<double>::Inf ();
                   value = 0;
                 }
               else
@@ -2992,7 +2992,7 @@
       // check for "treat as empty" string
       if (treat_as_empty.numel ()
           && (is.fail () || octave_is_NaN_or_NA (Complex (re))
-              || re == octave_Inf))
+              || re == octave_numeric_limits<double>::Inf ()))
         {
 
           for (int i = 0; i < treat_as_empty.numel (); i++)
--- a/libinterp/corefcn/quadcc.cc	Mon May 16 23:20:50 2016 -0400
+++ b/libinterp/corefcn/quadcc.cc	Tue May 17 01:11:27 2016 -0400
@@ -1702,7 +1702,7 @@
       Vinvfx (iv->fx, &(iv->c[idx[2]]), 2);
       Vinvfx (iv->fx, &(iv->c[0]), 0);
       for (i = 0; i < nnans; i++)
-        iv->fx[nans[i]] = octave_NaN;
+        iv->fx[nans[i]] = octave_numeric_limits<double>::NaN ();
       iv->a = iivals[j];
       iv->b = iivals[j + 1];
       iv->depth = 3;
@@ -1822,7 +1822,7 @@
             {
               downdate (&(iv->c[idx[d]]), n[d], d, nans, nnans);
               for (i = 0; i < nnans; i++)
-                iv->fx[nans[i]] = octave_NaN;
+                iv->fx[nans[i]] = octave_numeric_limits<double>::NaN ();
             }
 
           // Compute the error estimate.
@@ -1956,7 +1956,7 @@
             {
               downdate (ivl->c, n[0], 0, nans, nnans);
               for (i = 0; i < nnans; i++)
-                ivl->fx[nans[i]] = octave_NaN;
+                ivl->fx[nans[i]] = octave_numeric_limits<double>::NaN ();
             }
           for (i = 0; i <= n[d]; i++)
             {
@@ -1982,7 +1982,7 @@
                                   && ivl->c[0] / iv->c[0] > 2);
           if (ivl->ndiv > ndiv_max && 2 * ivl->ndiv > ivl->rdepth)
             {
-              igral = gnulib::copysign (octave_Inf, igral);
+              igral = gnulib::copysign (octave_numeric_limits<double>::Inf (), igral);
               warning ("quadcc: divergent integral detected");
               break;
             }
@@ -2047,7 +2047,7 @@
             {
               downdate (ivr->c, n[0], 0, nans, nnans);
               for (i = 0; i < nnans; i++)
-                ivr->fx[nans[i]] = octave_NaN;
+                ivr->fx[nans[i]] = octave_numeric_limits<double>::NaN ();
             }
           for (i = 0; i <= n[d]; i++)
             {
@@ -2073,7 +2073,7 @@
                                   && ivr->c[0] / iv->c[0] > 2);
           if (ivr->ndiv > ndiv_max && 2 * ivr->ndiv > ivr->rdepth)
             {
-              igral = gnulib::copysign (octave_Inf, igral);
+              igral = gnulib::copysign (octave_numeric_limits<double>::Inf (), igral);
               warning ("quadcc: divergent integral detected");
               break;
             }
--- a/libinterp/corefcn/sparse-xdiv.cc	Mon May 16 23:20:50 2016 -0400
+++ b/libinterp/corefcn/sparse-xdiv.cc	Tue May 17 01:11:27 2016 -0400
@@ -375,11 +375,11 @@
 
   Matrix result;
   if (a == 0.)
-    result = Matrix (nr, nc, octave_NaN);
+    result = Matrix (nr, nc, octave_numeric_limits<double>::NaN ());
   else if (a > 0.)
-    result = Matrix (nr, nc, octave_Inf);
+    result = Matrix (nr, nc, octave_numeric_limits<double>::Inf ());
   else
-    result = Matrix (nr, nc, -octave_Inf);
+    result = Matrix (nr, nc, -octave_numeric_limits<double>::Inf ());
 
 
   for (octave_idx_type j = 0; j < nc; j++)
@@ -398,7 +398,7 @@
   octave_idx_type nr = b.rows ();
   octave_idx_type nc = b.cols ();
 
-  ComplexMatrix result (nr, nc, Complex (octave_NaN, octave_NaN));
+  ComplexMatrix result (nr, nc, Complex (octave_numeric_limits<double>::NaN (), octave_numeric_limits<double>::NaN ()));
 
   for (octave_idx_type j = 0; j < nc; j++)
     for (octave_idx_type i = b.cidx (j); i < b.cidx (j+1); i++)
--- a/libinterp/corefcn/str2double.cc	Mon May 16 23:20:50 2016 -0400
+++ b/libinterp/corefcn/str2double.cc	Tue May 17 01:11:27 2016 -0400
@@ -64,7 +64,7 @@
       char c2 = is.get ();
       if (std::tolower (c1) == 'n' && std::tolower (c2) == 'f')
         {
-          num = octave_Inf;
+          num = octave_numeric_limits<double>::Inf ();
           is.peek (); // May set EOF bit.
         }
       else
@@ -85,7 +85,7 @@
           char c2 = is.get ();
           if (c1 == 'a' && c2 == 'N')
             {
-              num = octave_NaN;
+              num = octave_numeric_limits<double>::NaN ();
               is.peek (); // May set EOF bit.
             }
           else
@@ -275,9 +275,9 @@
   bool i1, i2, s1, s2;
 
   if (is.eof ())
-    val = octave_NaN;
+    val = octave_numeric_limits<double>::NaN ();
   else if (! extract_num (is, num, i1, s1))
-    val = octave_NaN;
+    val = octave_numeric_limits<double>::NaN ();
   else
     {
       set_component (val, num, i1);
@@ -285,7 +285,7 @@
       if (! is.eof ())
         {
           if (! extract_num (is, num, i2, s2) || i1 == i2 || ! s2)
-            val = octave_NaN;
+            val = octave_numeric_limits<double>::NaN ();
           else
             set_component (val, num, i2);
         }
@@ -349,7 +349,7 @@
   if (args(0).is_string ())
     {
       if (args(0).rows () == 0 || args(0).columns () == 0)
-        retval = Matrix (1, 1, octave_NaN);
+        retval = Matrix (1, 1, octave_numeric_limits<double>::NaN ());
       else if (args(0).rows () == 1 && args(0).ndims () == 2)
         retval = str2double1 (args(0).string_value ());
       else
@@ -363,7 +363,7 @@
     {
       const Cell cell = args(0).cell_value ();
 
-      ComplexNDArray output (cell.dims (), octave_NaN);
+      ComplexNDArray output (cell.dims (), octave_numeric_limits<double>::NaN ());
 
       for (octave_idx_type i = 0; i < cell.numel (); i++)
         {
@@ -373,7 +373,7 @@
       retval = output;
     }
   else
-    retval = Matrix (1, 1, octave_NaN);
+    retval = Matrix (1, 1, octave_numeric_limits<double>::NaN ());
 
   return retval;
 }
--- a/libinterp/corefcn/variables.h	Mon May 16 23:20:50 2016 -0400
+++ b/libinterp/corefcn/variables.h	Tue May 17 01:11:27 2016 -0400
@@ -108,8 +108,8 @@
 extern OCTINTERP_API octave_value
 set_internal_variable (double& var, const octave_value_list& args,
                        int nargout, const char *nm,
-                       double minval = -octave_Inf,
-                       double maxval = octave_Inf);
+                       double minval = -octave_numeric_limits<double>::Inf (),
+                       double maxval = octave_numeric_limits<double>::Inf ());
 
 extern OCTINTERP_API octave_value
 set_internal_variable (std::string& var, const octave_value_list& args,
--- a/libinterp/dldfcn/__glpk__.cc	Mon May 16 23:20:50 2016 -0400
+++ b/libinterp/dldfcn/__glpk__.cc	Tue May 17 01:11:27 2016 -0400
@@ -422,7 +422,7 @@
       if (xisinf (lb[i]))
         {
           freeLB(i) = 1;
-          lb[i] = -octave_Inf;
+          lb[i] = -octave_numeric_limits<double>::Inf ();
         }
       else
         freeLB(i) = 0;
@@ -443,7 +443,7 @@
       if (xisinf (ub[i]))
         {
           freeUB(i) = 1;
-          ub[i] = octave_Inf;
+          ub[i] = octave_numeric_limits<double>::Inf ();
         }
       else
         freeUB(i) = 0;
--- a/libinterp/dldfcn/__init_fltk__.cc	Mon May 16 23:20:50 2016 -0400
+++ b/libinterp/dldfcn/__init_fltk__.cc	Tue May 17 01:11:27 2016 -0400
@@ -1199,7 +1199,7 @@
                      int px1 = -1, int py1 = -1)
   {
     double x0, y0, x1, y1;
-    x0 = y0 = x1 = y1 = octave_NaN;
+    x0 = y0 = x1 = y1 = octave_numeric_limits<double>::NaN ();
     std::stringstream cbuf;
     cbuf.precision (4);
     cbuf.width (6);
--- a/libinterp/dldfcn/__voronoi__.cc	Mon May 16 23:20:50 2016 -0400
+++ b/libinterp/dldfcn/__voronoi__.cc	Tue May 17 01:11:27 2016 -0400
@@ -244,7 +244,7 @@
   Matrix F (num_voronoi_vertices+1, dim);
 
   for (octave_idx_type d = 0; d < dim; d++)
-    F(0,d) = octave_Inf;
+    F(0,d) = octave_numeric_limits<double>::Inf ();
 
   // The cell array of vectors of indices into F that represent the
   // vertices of the Voronoi regions (cells).
--- a/libinterp/octave.cc	Mon May 16 23:20:50 2016 -0400
+++ b/libinterp/octave.cc	Tue May 17 01:11:27 2016 -0400
@@ -809,9 +809,12 @@
 
   // The idea here is to force xerbla to be referenced so that we will link to
   // our own version instead of the one provided by the BLAS library.  But
-  // octave_NaN should never be -1, so we should never actually call xerbla.
+  // octave_numeric_limits<double>::NaN () should never be -1, so we
+  // should never actually call xerbla.  FIXME (again!):  If this
+  // becomes a constant expression the test might be optimized away and
+  // then the reference to the function might also disappear.
 
-  if (octave_NaN == -1)
+  if (octave_numeric_limits<double>::NaN () == -1)
     F77_FUNC (xerbla, XERBLA) ("octave", 13 F77_CHAR_ARG_LEN (6));
 
   initialize_error_handlers ();
--- a/liboctave/array/CMatrix.cc	Mon May 16 23:20:50 2016 -0400
+++ b/liboctave/array/CMatrix.cc	Tue May 17 01:11:27 2016 -0400
@@ -263,7 +263,7 @@
                                F77_CHAR_ARG_LEN_DECL);
 }
 
-static const Complex Complex_NaN_result (octave_NaN, octave_NaN);
+static const Complex Complex_NaN_result (octave_numeric_limits<double>::NaN (), octave_numeric_limits<double>::NaN ());
 
 // Complex Matrix class
 
@@ -1126,7 +1126,7 @@
         ret = finverse (mattype, info, rcon, force, calc_cond);
 
       if ((mattype.is_hermitian () || calc_cond) && rcon == 0.)
-        ret = ComplexMatrix (rows (), columns (), Complex (octave_Inf, 0.));
+        ret = ComplexMatrix (rows (), columns (), Complex (octave_numeric_limits<double>::Inf (), 0.));
     }
 
   return ret;
@@ -1675,7 +1675,7 @@
 double
 ComplexMatrix::rcond (MatrixType &mattype) const
 {
-  double rcon = octave_NaN;
+  double rcon = octave_numeric_limits<double>::NaN ();
   octave_idx_type nr = rows ();
   octave_idx_type nc = cols ();
 
@@ -1683,7 +1683,7 @@
     (*current_liboctave_error_handler) ("matrix must be square");
 
   if (nr == 0 || nc == 0)
-    rcon = octave_Inf;
+    rcon = octave_numeric_limits<double>::Inf ();
   else
     {
       volatile int typ = mattype.type ();
@@ -3147,7 +3147,7 @@
 
           Complex tmp_min;
 
-          double abs_min = octave_NaN;
+          double abs_min = octave_numeric_limits<double>::NaN ();
 
           for (idx_j = 0; idx_j < nc; idx_j++)
             {
@@ -3222,7 +3222,7 @@
 
           Complex tmp_max;
 
-          double abs_max = octave_NaN;
+          double abs_max = octave_numeric_limits<double>::NaN ();
 
           for (idx_j = 0; idx_j < nc; idx_j++)
             {
@@ -3297,7 +3297,7 @@
 
           Complex tmp_min;
 
-          double abs_min = octave_NaN;
+          double abs_min = octave_numeric_limits<double>::NaN ();
 
           for (idx_i = 0; idx_i < nr; idx_i++)
             {
@@ -3372,7 +3372,7 @@
 
           Complex tmp_max;
 
-          double abs_max = octave_NaN;
+          double abs_max = octave_numeric_limits<double>::NaN ();
 
           for (idx_i = 0; idx_i < nr; idx_i++)
             {
--- a/liboctave/array/CNDArray.cc	Mon May 16 23:20:50 2016 -0400
+++ b/liboctave/array/CNDArray.cc	Tue May 17 01:11:27 2016 -0400
@@ -668,7 +668,8 @@
   return retval;
 }
 
-static const Complex Complex_NaN_result (octave_NaN, octave_NaN);
+static const Complex Complex_NaN_result (octave_numeric_limits<double>::NaN (),
+                                         octave_numeric_limits<double>::NaN ());
 
 ComplexNDArray
 ComplexNDArray::max (int dim) const
--- a/liboctave/array/CSparse.cc	Mon May 16 23:20:50 2016 -0400
+++ b/liboctave/array/CSparse.cc	Tue May 17 01:11:27 2016 -0400
@@ -252,7 +252,7 @@
   return false;
 }
 
-static const Complex Complex_NaN_result (octave_NaN, octave_NaN);
+static const Complex Complex_NaN_result (octave_numeric_limits<double>::NaN (), octave_numeric_limits<double>::NaN ());
 
 SparseComplexMatrix
 SparseComplexMatrix::max (int dim) const
@@ -290,7 +290,7 @@
       for (octave_idx_type j = 0; j < nc; j++)
         {
           Complex tmp_max;
-          double abs_max = octave_NaN;
+          double abs_max = octave_numeric_limits<double>::NaN ();
           octave_idx_type idx_j = 0;
           for (octave_idx_type i = cidx (j); i < cidx (j+1); i++)
             {
@@ -447,7 +447,7 @@
       for (octave_idx_type j = 0; j < nc; j++)
         {
           Complex tmp_min;
-          double abs_min = octave_NaN;
+          double abs_min = octave_numeric_limits<double>::NaN ();
           octave_idx_type idx_j = 0;
           for (octave_idx_type i = cidx (j); i < cidx (j+1); i++)
             {
@@ -785,7 +785,7 @@
   if (calccond)
     {
       double dmax = 0.;
-      double dmin = octave_Inf;
+      double dmin = octave_numeric_limits<double>::Inf ();
       for (octave_idx_type i = 0; i < nr; i++)
         {
           double tmp = std::abs (v[i]);
@@ -1308,7 +1308,7 @@
       if (calc_cond)
         {
           double dmax = 0.;
-          double dmin = octave_Inf;
+          double dmin = octave_numeric_limits<double>::Inf ();
           for (octave_idx_type i = 0; i < nm; i++)
             {
               double tmp = std::abs (data (i));
@@ -1398,7 +1398,7 @@
       if (calc_cond)
         {
           double dmax = 0.;
-          double dmin = octave_Inf;
+          double dmin = octave_numeric_limits<double>::Inf ();
           for (octave_idx_type i = 0; i < nm; i++)
             {
               double tmp = std::abs (data (i));
@@ -1458,7 +1458,7 @@
       if (calc_cond)
         {
           double dmax = 0.;
-          double dmin = octave_Inf;
+          double dmin = octave_numeric_limits<double>::Inf ();
           for (octave_idx_type i = 0; i < nr; i++)
             {
               double tmp = std::abs (data (i));
@@ -1548,7 +1548,7 @@
       if (calc_cond)
         {
           double dmax = 0.;
-          double dmin = octave_Inf;
+          double dmin = octave_numeric_limits<double>::Inf ();
           for (octave_idx_type i = 0; i < nm; i++)
             {
               double tmp = std::abs (data (i));
--- a/liboctave/array/Sparse-C.cc	Mon May 16 23:20:50 2016 -0400
+++ b/liboctave/array/Sparse-C.cc	Tue May 17 01:11:27 2016 -0400
@@ -37,7 +37,8 @@
 static double
 xabs (const Complex& x)
 {
-  return (xisinf (x.real ()) || xisinf (x.imag ())) ? octave_Inf : abs (x);
+  return ((xisinf (x.real ()) || xisinf (x.imag ()))
+          ? octave_numeric_limits<double>::Inf () : abs (x));
 }
 
 
--- a/liboctave/array/dDiagMatrix.cc	Mon May 16 23:20:50 2016 -0400
+++ b/liboctave/array/dDiagMatrix.cc	Tue May 17 01:11:27 2016 -0400
@@ -244,7 +244,7 @@
   for (octave_idx_type i = 0; i < len; i++)
     {
       if (elem (i, i) == 0.0)
-        retval.elem (i, i) = octave_Inf;
+        retval.elem (i, i) = octave_numeric_limits<double>::Inf ();
       else
         retval.elem (i, i) = 1.0 / elem (i, i);
     }
--- a/liboctave/array/dMatrix.cc	Mon May 16 23:20:50 2016 -0400
+++ b/liboctave/array/dMatrix.cc	Tue May 17 01:11:27 2016 -0400
@@ -816,7 +816,7 @@
         ret = finverse (mattype, info, rcon, force, calc_cond);
 
       if ((mattype.is_hermitian () || calc_cond) && rcon == 0.)
-        ret = Matrix (rows (), columns (), octave_Inf);
+        ret = Matrix (rows (), columns (), octave_numeric_limits<double>::Inf ());
     }
 
   return ret;
@@ -1352,7 +1352,7 @@
 double
 Matrix::rcond (MatrixType &mattype) const
 {
-  double rcon = octave_NaN;
+  double rcon = octave_numeric_limits<double>::NaN ();
   octave_idx_type nr = rows ();
   octave_idx_type nc = cols ();
 
@@ -1360,7 +1360,7 @@
     (*current_liboctave_error_handler) ("matrix must be square");
 
   if (nr == 0 || nc == 0)
-    rcon = octave_Inf;
+    rcon = octave_numeric_limits<double>::Inf ();
   else
     {
       volatile int typ = mattype.type ();
@@ -2693,7 +2693,7 @@
         {
           octave_idx_type idx_j;
 
-          double tmp_min = octave_NaN;
+          double tmp_min = octave_numeric_limits<double>::NaN ();
 
           for (idx_j = 0; idx_j < nc; idx_j++)
             {
@@ -2748,7 +2748,7 @@
         {
           octave_idx_type idx_j;
 
-          double tmp_max = octave_NaN;
+          double tmp_max = octave_numeric_limits<double>::NaN ();
 
           for (idx_j = 0; idx_j < nc; idx_j++)
             {
@@ -2803,7 +2803,7 @@
         {
           octave_idx_type idx_i;
 
-          double tmp_min = octave_NaN;
+          double tmp_min = octave_numeric_limits<double>::NaN ();
 
           for (idx_i = 0; idx_i < nr; idx_i++)
             {
@@ -2858,7 +2858,7 @@
         {
           octave_idx_type idx_i;
 
-          double tmp_max = octave_NaN;
+          double tmp_max = octave_numeric_limits<double>::NaN ();
 
           for (idx_i = 0; idx_i < nr; idx_i++)
             {
--- a/liboctave/array/dSparse.cc	Mon May 16 23:20:50 2016 -0400
+++ b/liboctave/array/dSparse.cc	Tue May 17 01:11:27 2016 -0400
@@ -302,7 +302,7 @@
       octave_idx_type nel = 0;
       for (octave_idx_type j = 0; j < nc; j++)
         {
-          double tmp_max = octave_NaN;
+          double tmp_max = octave_numeric_limits<double>::NaN ();
           octave_idx_type idx_j = 0;
           for (octave_idx_type i = cidx (j); i < cidx (j+1); i++)
             {
@@ -401,7 +401,7 @@
           if (idx_arg(j) == -1)
             {
               idx_arg(j) = 0;
-              result.xdata (ii) = octave_NaN;
+              result.xdata (ii) = octave_numeric_limits<double>::NaN ();
               result.xridx (ii++) = j;
             }
           else
@@ -453,7 +453,7 @@
       octave_idx_type nel = 0;
       for (octave_idx_type j = 0; j < nc; j++)
         {
-          double tmp_min = octave_NaN;
+          double tmp_min = octave_numeric_limits<double>::NaN ();
           octave_idx_type idx_j = 0;
           for (octave_idx_type i = cidx (j); i < cidx (j+1); i++)
             {
@@ -552,7 +552,7 @@
           if (idx_arg(j) == -1)
             {
               idx_arg(j) = 0;
-              result.xdata (ii) = octave_NaN;
+              result.xdata (ii) = octave_numeric_limits<double>::NaN ();
               result.xridx (ii++) = j;
             }
           else
@@ -876,7 +876,7 @@
   if (calccond)
     {
       double dmax = 0.;
-      double dmin = octave_Inf;
+      double dmin = octave_numeric_limits<double>::Inf ();
       for (octave_idx_type i = 0; i < nr; i++)
         {
           double tmp = fabs (v[i]);
@@ -1390,7 +1390,7 @@
       if (calc_cond)
         {
           double dmax = 0.;
-          double dmin = octave_Inf;
+          double dmin = octave_numeric_limits<double>::Inf ();
           for (octave_idx_type i = 0; i < nm; i++)
             {
               double tmp = fabs (data (i));
@@ -1479,7 +1479,7 @@
       if (calc_cond)
         {
           double dmax = 0.;
-          double dmin = octave_Inf;
+          double dmin = octave_numeric_limits<double>::Inf ();
           for (octave_idx_type i = 0; i < nm; i++)
             {
               double tmp = fabs (data (i));
@@ -1538,7 +1538,7 @@
       if (calc_cond)
         {
           double dmax = 0.;
-          double dmin = octave_Inf;
+          double dmin = octave_numeric_limits<double>::Inf ();
           for (octave_idx_type i = 0; i < nm; i++)
             {
               double tmp = fabs (data (i));
@@ -1627,7 +1627,7 @@
       if (calc_cond)
         {
           double dmax = 0.;
-          double dmin = octave_Inf;
+          double dmin = octave_numeric_limits<double>::Inf ();
           for (octave_idx_type i = 0; i < nm; i++)
             {
               double tmp = fabs (data (i));
--- a/liboctave/array/fCMatrix.cc	Mon May 16 23:20:50 2016 -0400
+++ b/liboctave/array/fCMatrix.cc	Tue May 17 01:11:27 2016 -0400
@@ -263,8 +263,8 @@
                                F77_CHAR_ARG_LEN_DECL);
 }
 
-static const FloatComplex FloatComplex_NaN_result (octave_Float_NaN,
-                                                   octave_Float_NaN);
+static const FloatComplex FloatComplex_NaN_result (octave_numeric_limits<float>::NaN (),
+                                                   octave_numeric_limits<float>::NaN ());
 
 // FloatComplex Matrix class
 
@@ -1130,7 +1130,7 @@
 
       if ((mattype.is_hermitian () || calc_cond) && rcon == 0.)
         ret = FloatComplexMatrix (rows (), columns (),
-                                  FloatComplex (octave_Float_Inf, 0.));
+                                  FloatComplex (octave_numeric_limits<float>::Inf (), 0.));
     }
 
   return ret;
@@ -1671,7 +1671,7 @@
 float
 FloatComplexMatrix::rcond (MatrixType &mattype) const
 {
-  float rcon = octave_NaN;
+  float rcon = octave_numeric_limits<float>::NaN ();
   octave_idx_type nr = rows ();
   octave_idx_type nc = cols ();
 
@@ -1679,7 +1679,7 @@
     (*current_liboctave_error_handler) ("matrix must be square");
 
   if (nr == 0 || nc == 0)
-    rcon = octave_Inf;
+    rcon = octave_numeric_limits<float>::Inf ();
   else
     {
       volatile int typ = mattype.type ();
@@ -3160,7 +3160,7 @@
 
           FloatComplex tmp_min;
 
-          float abs_min = octave_Float_NaN;
+          float abs_min = octave_numeric_limits<float>::NaN ();
 
           for (idx_j = 0; idx_j < nc; idx_j++)
             {
@@ -3235,7 +3235,7 @@
 
           FloatComplex tmp_max;
 
-          float abs_max = octave_Float_NaN;
+          float abs_max = octave_numeric_limits<float>::NaN ();
 
           for (idx_j = 0; idx_j < nc; idx_j++)
             {
@@ -3310,7 +3310,7 @@
 
           FloatComplex tmp_min;
 
-          float abs_min = octave_Float_NaN;
+          float abs_min = octave_numeric_limits<float>::NaN ();
 
           for (idx_i = 0; idx_i < nr; idx_i++)
             {
@@ -3385,7 +3385,7 @@
 
           FloatComplex tmp_max;
 
-          float abs_max = octave_Float_NaN;
+          float abs_max = octave_numeric_limits<float>::NaN ();
 
           for (idx_i = 0; idx_i < nr; idx_i++)
             {
--- a/liboctave/array/fCNDArray.cc	Mon May 16 23:20:50 2016 -0400
+++ b/liboctave/array/fCNDArray.cc	Tue May 17 01:11:27 2016 -0400
@@ -675,8 +675,8 @@
   return retval;
 }
 
-static const FloatComplex FloatComplex_NaN_result (octave_Float_NaN,
-                                                   octave_Float_NaN);
+static const FloatComplex FloatComplex_NaN_result (octave_numeric_limits<float>::NaN (),
+                                                   octave_numeric_limits<float>::NaN ());
 
 FloatComplexNDArray
 FloatComplexNDArray::max (int dim) const
--- a/liboctave/array/fDiagMatrix.cc	Mon May 16 23:20:50 2016 -0400
+++ b/liboctave/array/fDiagMatrix.cc	Tue May 17 01:11:27 2016 -0400
@@ -244,7 +244,7 @@
   for (octave_idx_type i = 0; i < len; i++)
     {
       if (elem (i, i) == 0.0)
-        retval.elem (i, i) = octave_Inf;
+        retval.elem (i, i) = octave_numeric_limits<float>::Inf ();
       else
         retval.elem (i, i) = 1.0 / elem (i, i);
     }
--- a/liboctave/array/fMatrix.cc	Mon May 16 23:20:50 2016 -0400
+++ b/liboctave/array/fMatrix.cc	Tue May 17 01:11:27 2016 -0400
@@ -823,7 +823,7 @@
         ret = finverse (mattype, info, rcon, force, calc_cond);
 
       if ((mattype.is_hermitian () || calc_cond) && rcon == 0.)
-        ret = FloatMatrix (rows (), columns (), octave_Float_Inf);
+        ret = FloatMatrix (rows (), columns (), octave_numeric_limits<float>::Inf ());
     }
 
   return ret;
@@ -1363,7 +1363,7 @@
 float
 FloatMatrix::rcond (MatrixType &mattype) const
 {
-  float rcon = octave_NaN;
+  float rcon = octave_numeric_limits<float>::NaN ();
   octave_idx_type nr = rows ();
   octave_idx_type nc = cols ();
 
@@ -1371,7 +1371,7 @@
     (*current_liboctave_error_handler) ("matrix must be square");
 
   if (nr == 0 || nc == 0)
-    rcon = octave_Inf;
+    rcon = octave_numeric_limits<float>::Inf ();
   else
     {
       volatile int typ = mattype.type ();
@@ -2710,7 +2710,7 @@
         {
           octave_idx_type idx_j;
 
-          float tmp_min = octave_Float_NaN;
+          float tmp_min = octave_numeric_limits<float>::NaN ();
 
           for (idx_j = 0; idx_j < nc; idx_j++)
             {
@@ -2765,7 +2765,7 @@
         {
           octave_idx_type idx_j;
 
-          float tmp_max = octave_Float_NaN;
+          float tmp_max = octave_numeric_limits<float>::NaN ();
 
           for (idx_j = 0; idx_j < nc; idx_j++)
             {
@@ -2820,7 +2820,7 @@
         {
           octave_idx_type idx_i;
 
-          float tmp_min = octave_Float_NaN;
+          float tmp_min = octave_numeric_limits<float>::NaN ();
 
           for (idx_i = 0; idx_i < nr; idx_i++)
             {
@@ -2875,7 +2875,7 @@
         {
           octave_idx_type idx_i;
 
-          float tmp_max = octave_Float_NaN;
+          float tmp_max = octave_numeric_limits<float>::NaN ();
 
           for (idx_i = 0; idx_i < nr; idx_i++)
             {
--- a/liboctave/cruft/Faddeeva/Faddeeva.cc	Mon May 16 23:20:50 2016 -0400
+++ b/liboctave/cruft/Faddeeva/Faddeeva.cc	Tue May 17 01:11:27 2016 -0400
@@ -166,8 +166,8 @@
 #  include <limits>
 
 // use std::numeric_limits, since 1./0. and 0./0. fail with some compilers (MS)
-#  define Inf octave_Inf
-#  define NaN octave_NaN
+#  define Inf octave_numeric_limits<double>::Inf ()
+#  define NaN octave_numeric_limits<double>::NaN ()
 
 typedef std::complex<double> cmplx;
 
--- a/liboctave/numeric/eigs-base.cc	Mon May 16 23:20:50 2016 -0400
+++ b/liboctave/numeric/eigs-base.cc	Tue May 17 01:11:27 2016 -0400
@@ -441,8 +441,8 @@
     }
 
   // Test condition number of LU decomposition
-  double minU = octave_NaN;
-  double maxU = octave_NaN;
+  double minU = octave_numeric_limits<double>::NaN ();
+  double maxU = octave_numeric_limits<double>::NaN ();
   for (octave_idx_type j = 0; j < n; j++)
     {
       double d = 0.;
@@ -517,8 +517,8 @@
     P[j] = Q[j] = j;
 
   // Test condition number of LU decomposition
-  double minU = octave_NaN;
-  double maxU = octave_NaN;
+  double minU = octave_numeric_limits<double>::NaN ();
+  double maxU = octave_numeric_limits<double>::NaN ();
   for (octave_idx_type j = 0; j < n; j++)
     {
       double d = std::abs (U.xelem (j,j));
@@ -605,8 +605,8 @@
     }
 
   // Test condition number of LU decomposition
-  double minU = octave_NaN;
-  double maxU = octave_NaN;
+  double minU = octave_numeric_limits<double>::NaN ();
+  double maxU = octave_numeric_limits<double>::NaN ();
   for (octave_idx_type j = 0; j < n; j++)
     {
       double d = 0.;
@@ -681,8 +681,8 @@
     P[j] = Q[j] = j;
 
   // Test condition number of LU decomposition
-  double minU = octave_NaN;
-  double maxU = octave_NaN;
+  double minU = octave_numeric_limits<double>::NaN ();
+  double maxU = octave_numeric_limits<double>::NaN ();
   for (octave_idx_type j = 0; j < n; j++)
     {
       double d = std::abs (U.xelem (j,j));
@@ -864,7 +864,7 @@
           // a value in this array to NaN and testing for it
           // is a way of obtaining the iteration counter.
           if (ido != 99)
-            workl[iptr(5)-1] = octave_NaN;
+            workl[iptr(5)-1] = octave_numeric_limits<double>::NaN ();
         }
 
       if (ido == -1 || ido == 1 || ido == 2)
@@ -1117,7 +1117,7 @@
           // a value in this array to NaN and testing for it
           // is a way of obtaining the iteration counter.
           if (ido != 99)
-            workl[iptr(5)-1] = octave_NaN;
+            workl[iptr(5)-1] = octave_numeric_limits<double>::NaN ();
         }
 
       if (ido == -1 || ido == 1 || ido == 2)
@@ -1390,7 +1390,7 @@
           // a value in this array to NaN and testing for it
           // is a way of obtaining the iteration counter.
           if (ido != 99)
-            workl[iptr(5)-1] = octave_NaN;
+            workl[iptr(5)-1] = octave_numeric_limits<double>::NaN ();
         }
 
 
@@ -1659,7 +1659,7 @@
           // a value in this array to NaN and testing for it
           // is a way of obtaining the iteration counter.
           if (ido != 99)
-            workl[iptr(5)-1] = octave_NaN;
+            workl[iptr(5)-1] = octave_numeric_limits<double>::NaN ();
         }
 
       if (ido == -1 || ido == 1 || ido == 2)
@@ -1961,7 +1961,7 @@
           // a value in this array to NaN and testing for it
           // is a way of obtaining the iteration counter.
           if (ido != 99)
-            workl[iptr(5)-1] = octave_NaN;
+            workl[iptr(5)-1] = octave_numeric_limits<double>::NaN ();
         }
 
       if (ido == -1 || ido == 1 || ido == 2)
@@ -2289,7 +2289,7 @@
           // a value in this array to NaN and testing for it
           // is a way of obtaining the iteration counter.
           if (ido != 99)
-            workl[iptr(5)-1] = octave_NaN;
+            workl[iptr(5)-1] = octave_numeric_limits<double>::NaN ();
         }
 
       if (ido == -1 || ido == 1 || ido == 2)
@@ -2609,7 +2609,7 @@
           // a value in this array to NaN and testing for it
           // is a way of obtaining the iteration counter.
           if (ido != 99)
-            workl[iptr(5)-1] = octave_NaN;
+            workl[iptr(5)-1] = octave_numeric_limits<double>::NaN ();
         }
 
       if (ido == -1 || ido == 1 || ido == 2)
@@ -2865,7 +2865,7 @@
           // a value in this array to NaN and testing for it
           // is a way of obtaining the iteration counter.
           if (ido != 99)
-            workl[iptr(5)-1] = octave_NaN;
+            workl[iptr(5)-1] = octave_numeric_limits<double>::NaN ();
         }
 
       if (ido == -1 || ido == 1 || ido == 2)
@@ -3150,7 +3150,7 @@
           // a value in this array to NaN and testing for it
           // is a way of obtaining the iteration counter.
           if (ido != 99)
-            workl[iptr(5)-1] = octave_NaN;
+            workl[iptr(5)-1] = octave_numeric_limits<double>::NaN ();
         }
 
       if (ido == -1 || ido == 1 || ido == 2)
--- a/liboctave/numeric/lo-mappers.cc	Mon May 16 23:20:50 2016 -0400
+++ b/liboctave/numeric/lo-mappers.cc	Tue May 17 01:11:27 2016 -0400
@@ -86,7 +86,7 @@
   else if (x > 0.0)
     tmp = 1.0;
 
-  return xisnan (x) ? octave_NaN : tmp;
+  return xisnan (x) ? octave_numeric_limits<double>::NaN () : tmp;
 }
 
 double
@@ -309,7 +309,7 @@
   else if (x > 0.0)
     tmp = 1.0;
 
-  return xisnan (x) ? octave_Float_NaN : tmp;
+  return xisnan (x) ? octave_numeric_limits<float>::NaN () : tmp;
 }
 
 float
--- a/liboctave/numeric/lo-mappers.h	Mon May 16 23:20:50 2016 -0400
+++ b/liboctave/numeric/lo-mappers.h	Tue May 17 01:11:27 2016 -0400
@@ -314,7 +314,7 @@
   T retval;
 
   if (y == 0)
-    retval = octave_NaN;
+    retval = octave_numeric_limits<T>::NaN ();
   else
     {
       T q = x / y;
--- a/liboctave/numeric/lo-specfun.cc	Mon May 16 23:20:50 2016 -0400
+++ b/liboctave/numeric/lo-specfun.cc	Tue May 17 01:11:27 2016 -0400
@@ -441,11 +441,13 @@
   // tgamma.  Matlab does not have -0.
 
   if (x == 0)
-    result = xnegative_sign (x) ? -octave_Inf : octave_Inf;
+    result = (xnegative_sign (x)
+              ? -octave_numeric_limits<double>::Inf ()
+              : octave_numeric_limits<double>::Inf ());
   else if ((x < 0 && D_NINT (x) == x) || xisinf (x))
-    result = octave_Inf;
+    result = octave_numeric_limits<double>::Inf ();
   else if (xisnan (x))
-    result = octave_NaN;
+    result = octave_numeric_limits<double>::NaN ();
   else
     {
 #if defined (HAVE_TGAMMA)
@@ -470,7 +472,7 @@
   if (xisnan (x))
     result = x;
   else if ((x <= 0 && D_NINT (x) == x) || xisinf (x))
-    result = octave_Inf;
+    result = octave_numeric_limits<double>::Inf ();
   else
     F77_XFCN (dlgams, DLGAMS, (x, result, sgngam));
 
@@ -492,7 +494,7 @@
   if (xisnan (x))
     result = x;
   else if ((x <= 0 && D_NINT (x) == x) || xisinf (x))
-    result = octave_Inf;
+    result = octave_numeric_limits<double>::Inf ();
   else
     F77_XFCN (dlgams, DLGAMS, (x, result, sgngam));
 
@@ -513,11 +515,13 @@
   // tgamma.  Matlab does not have -0.
 
   if (x == 0)
-    result = xnegative_sign (x) ? -octave_Float_Inf : octave_Float_Inf;
+    result = (xnegative_sign (x)
+              ? -octave_numeric_limits<float>::Inf ()
+              : octave_numeric_limits<float>::Inf ());
   else if ((x < 0 && D_NINT (x) == x) || xisinf (x))
-    result = octave_Float_Inf;
+    result = octave_numeric_limits<float>::Inf ();
   else if (xisnan (x))
-    result = octave_Float_NaN;
+    result = octave_numeric_limits<float>::NaN ();
   else
     {
 #if defined (HAVE_TGAMMA)
@@ -542,7 +546,7 @@
   if (xisnan (x))
     result = x;
   else if ((x <= 0 && D_NINT (x) == x) || xisinf (x))
-    result = octave_Float_Inf;
+    result = octave_numeric_limits<float>::Inf ();
   else
     F77_XFCN (algams, ALGAMS, (x, result, sgngam));
 
@@ -564,7 +568,7 @@
   if (xisnan (x))
     result = x;
   else if ((x <= 0 && D_NINT (x) == x) || xisinf (x))
-    result = octave_Float_Inf;
+    result = octave_numeric_limits<float>::Inf ();
   else
     F77_XFCN (algams, ALGAMS, (x, result, sgngam));
 
@@ -837,8 +841,13 @@
 static inline Complex
 bessel_return_value (const Complex& val, octave_idx_type ierr)
 {
-  static const Complex inf_val = Complex (octave_Inf, octave_Inf);
-  static const Complex nan_val = Complex (octave_NaN, octave_NaN);
+  static const Complex inf_val
+    = Complex (octave_numeric_limits<double>::Inf (),
+               octave_numeric_limits<double>::Inf ());
+
+  static const Complex nan_val
+    = Complex (octave_numeric_limits<double>::NaN (),
+               octave_numeric_limits<double>::NaN ());
 
   Complex retval;
 
@@ -918,7 +927,8 @@
           retval = bessel_return_value (tmp, ierr);
         }
       else
-        retval = Complex (octave_NaN, octave_NaN);
+        retval = Complex (octave_numeric_limits<double>::NaN (),
+                          octave_numeric_limits<double>::NaN ());
     }
 
   return retval;
@@ -945,7 +955,7 @@
 
       if (zr == 0.0 && zi == 0.0)
         {
-          yr = -octave_Inf;
+          yr = -octave_numeric_limits<double>::Inf ();
           yi = 0.0;
         }
       else
@@ -988,7 +998,8 @@
           retval = bessel_return_value (tmp, ierr);
         }
       else
-        retval = Complex (octave_NaN, octave_NaN);
+        retval = Complex (octave_numeric_limits<double>::NaN (),
+                          octave_numeric_limits<double>::NaN ());
     }
 
   return retval;
@@ -1052,7 +1063,8 @@
           retval = bessel_return_value (tmp, ierr);
         }
       else
-        retval = Complex (octave_NaN, octave_NaN);
+        retval = Complex (octave_numeric_limits<double>::NaN (),
+                          octave_numeric_limits<double>::NaN ());
     }
 
   return retval;
@@ -1077,7 +1089,7 @@
 
       if (zr == 0.0 && zi == 0.0)
         {
-          yr = octave_Inf;
+          yr = octave_numeric_limits<double>::Inf ();
           yi = 0.0;
         }
       else
@@ -1468,10 +1480,13 @@
 static inline FloatComplex
 bessel_return_value (const FloatComplex& val, octave_idx_type ierr)
 {
-  static const FloatComplex inf_val = FloatComplex (octave_Float_Inf,
-                                                    octave_Float_Inf);
-  static const FloatComplex nan_val = FloatComplex (octave_Float_NaN,
-                                                    octave_Float_NaN);
+  static const FloatComplex inf_val
+    = FloatComplex (octave_numeric_limits<float>::Inf (),
+                    octave_numeric_limits<float>::Inf ());
+
+  static const FloatComplex nan_val
+    = FloatComplex (octave_numeric_limits<float>::NaN (),
+                    octave_numeric_limits<float>::NaN ());
 
   FloatComplex retval;
 
@@ -1548,7 +1563,8 @@
           retval = bessel_return_value (tmp, ierr);
         }
       else
-        retval = FloatComplex (octave_Float_NaN, octave_Float_NaN);
+        retval = FloatComplex (octave_numeric_limits<float>::NaN (),
+                               octave_numeric_limits<float>::NaN ());
     }
 
   return retval;
@@ -1571,7 +1587,7 @@
 
       if (real (z) == 0.0 && imag (z) == 0.0)
         {
-          y = FloatComplex (-octave_Float_Inf, 0.0);
+          y = FloatComplex (-octave_numeric_limits<float>::Inf (), 0.0);
         }
       else
         {
@@ -1613,7 +1629,8 @@
           retval = bessel_return_value (tmp, ierr);
         }
       else
-        retval = FloatComplex (octave_Float_NaN, octave_Float_NaN);
+        retval = FloatComplex (octave_numeric_limits<float>::NaN (),
+                               octave_numeric_limits<float>::NaN ());
     }
 
   return retval;
@@ -1666,7 +1683,8 @@
           retval = bessel_return_value (tmp, ierr);
         }
       else
-        retval = FloatComplex (octave_Float_NaN, octave_Float_NaN);
+        retval = FloatComplex (octave_numeric_limits<float>::NaN (),
+                               octave_numeric_limits<float>::NaN ());
     }
 
   return retval;
@@ -1687,7 +1705,7 @@
 
       if (real (z) == 0.0 && imag (z) == 0.0)
         {
-          y = FloatComplex (octave_Float_Inf, 0.0);
+          y = FloatComplex (octave_numeric_limits<float>::Inf (), 0.0);
         }
       else
         {
@@ -3031,9 +3049,9 @@
       y = yn / yd * signum (-x);
     }
   else if (ax == 1.0)
-    return octave_Inf * signum (x);
+    return octave_numeric_limits<double>::Inf () * signum (x);
   else
-    return octave_NaN;
+    return octave_numeric_limits<double>::NaN ();
 
   if (refine)
     {
@@ -3118,11 +3136,11 @@
         y = -y;
     }
   else if (x == 0.0)
-    return octave_Inf;
+    return octave_numeric_limits<double>::Inf ();
   else if (x == 2.0)
-    return -octave_Inf;
+    return -octave_numeric_limits<double>::Inf ();
   else
-    return octave_NaN;
+    return octave_numeric_limits<double>::NaN ();
 
   if (refine)
     {
@@ -3756,7 +3774,7 @@
     {
       // limits - zeros of the gamma function
       if (is_int)
-        p = -octave_Inf; // Matlab returns -Inf for psi (0)
+        p = -octave_numeric_limits<T>::Inf (); // Matlab returns -Inf for psi (0)
       else
         // Abramowitz and Stegun, page 259, eq 6.3.7
         p = psi (1 - z) - (pi / tan (pi * z));
@@ -3892,9 +3910,9 @@
         ans = -ans;
     }
   else if (ierr == 2)
-    ans = - octave_Inf;
+    ans = - octave_numeric_limits<T>::Inf ();
   else // we probably never get here
-    ans = octave_NaN;
+    ans = octave_numeric_limits<T>::NaN ();
 
   return ans;
 }
--- a/liboctave/numeric/oct-norm.cc	Mon May 16 23:20:50 2016 -0400
+++ b/liboctave/numeric/oct-norm.cc	Tue May 17 01:11:27 2016 -0400
@@ -180,7 +180,7 @@
   void accum (U val)
   {
     if (xisnan (val))
-      max = octave_NaN;
+      max = octave_numeric_limits<R>::NaN ();
     else
       max = std::max (max, std::abs (val));
   }
@@ -193,12 +193,12 @@
 {
   R min;
 public:
-  norm_accumulator_minf () : min (octave_Inf) {}
+  norm_accumulator_minf () : min (octave_numeric_limits<R>::Inf ()) {}
   template <typename U>
   void accum (U val)
   {
     if (xisnan (val))
-      min = octave_NaN;
+      min = octave_numeric_limits<R>::NaN ();
     else
       min = std::min (min, std::abs (val));
   }
--- a/liboctave/numeric/oct-rand.cc	Mon May 16 23:20:50 2016 -0400
+++ b/liboctave/numeric/oct-rand.cc	Tue May 17 01:11:27 2016 -0400
@@ -357,7 +357,7 @@
 
         case poisson_dist:
           if (a < 0.0 || ! xfinite (a))
-            retval = octave_NaN;
+            retval = octave_numeric_limits<double>::NaN ();
           else
             {
               // workaround bug in ignpoi, by calling with different Mu
@@ -368,7 +368,7 @@
 
         case gamma_dist:
           if (a <= 0.0 || ! xfinite (a))
-            retval = octave_NaN;
+            retval = octave_numeric_limits<double>::NaN ();
           else
             F77_FUNC (dgengam, DGENGAM) (1.0, a, retval);
           break;
@@ -440,7 +440,7 @@
 
         case poisson_dist:
           if (da < 0.0 || ! xfinite (a))
-            dretval = octave_NaN;
+            dretval = octave_numeric_limits<double>::NaN ();
           else
             {
               // workaround bug in ignpoi, by calling with different Mu
@@ -451,7 +451,7 @@
 
         case gamma_dist:
           if (da <= 0.0 || ! xfinite (a))
-            dretval = octave_NaN;
+            dretval = octave_numeric_limits<double>::NaN ();
           else
             F77_FUNC (dgengam, DGENGAM) (1.0, da, dretval);
           break;
@@ -769,7 +769,7 @@
       if (use_old_generators)
         {
           if (a < 0.0 || ! xfinite (a))
-#define RAND_FUNC(x) x = octave_NaN;
+#define RAND_FUNC(x) x = octave_numeric_limits<double>::NaN ();
             MAKE_RAND (len);
 #undef RAND_FUNC
           else
@@ -790,7 +790,7 @@
       if (use_old_generators)
         {
           if (a <= 0.0 || ! xfinite (a))
-#define RAND_FUNC(x) x = octave_NaN;
+#define RAND_FUNC(x) x = octave_numeric_limits<double>::NaN ();
             MAKE_RAND (len);
 #undef RAND_FUNC
           else
@@ -859,7 +859,7 @@
         {
           double da = a;
           if (da < 0.0 || ! xfinite (a))
-#define RAND_FUNC(x) x = octave_NaN;
+#define RAND_FUNC(x) x = octave_numeric_limits<double>::NaN ();
             MAKE_RAND (len);
 #undef RAND_FUNC
           else
@@ -881,7 +881,7 @@
         {
           double da = a;
           if (da <= 0.0 || ! xfinite (a))
-#define RAND_FUNC(x) x = octave_NaN;
+#define RAND_FUNC(x) x = octave_numeric_limits<double>::NaN ();
             MAKE_RAND (len);
 #undef RAND_FUNC
           else
--- a/liboctave/numeric/oct-rand.h	Mon May 16 23:20:50 2016 -0400
+++ b/liboctave/numeric/oct-rand.h	Tue May 17 01:11:27 2016 -0400
@@ -50,7 +50,7 @@
   // Return the current seed.
   static double seed (void)
   {
-    return instance_ok () ? instance->do_seed () : octave_NaN;
+    return instance_ok () ? instance->do_seed () : octave_numeric_limits<double>::NaN ();
   }
 
   // Set the seed.
@@ -135,13 +135,13 @@
   // Return the next number from the sequence.
   static double scalar (double a = 1.0)
   {
-    return instance_ok () ? instance->do_scalar (a) : octave_NaN;
+    return instance_ok () ? instance->do_scalar (a) : octave_numeric_limits<double>::NaN ();
   }
 
   // Return the next number from the sequence.
   static float float_scalar (float a = 1.0)
   {
-    return instance_ok () ? instance->do_float_scalar (a) : octave_Float_NaN;
+    return instance_ok () ? instance->do_float_scalar (a) : octave_numeric_limits<float>::NaN ();
   }
 
   // Return an array of numbers from the sequence.
--- a/liboctave/numeric/oct-spparms.cc	Mon May 16 23:20:50 2016 -0400
+++ b/liboctave/numeric/oct-spparms.cc	Tue May 17 01:11:27 2016 -0400
@@ -94,7 +94,7 @@
 double
 octave_sparse_params::get_key (const std::string& key)
 {
-  return instance_ok () ? instance->do_get_key (key) : octave_NaN;
+  return instance_ok () ? instance->do_get_key (key) : octave_numeric_limits<double>::NaN ();
 }
 
 double
@@ -209,7 +209,7 @@
         return params(i);
     }
 
-  return octave_NaN;
+  return octave_numeric_limits<double>::NaN ();
 }
 
 void
--- a/liboctave/numeric/randgamma.cc	Mon May 16 23:20:50 2016 -0400
+++ b/liboctave/numeric/randgamma.cc	Tue May 17 01:11:27 2016 -0400
@@ -86,8 +86,6 @@
 #include "randmtzig.h"
 #include "randgamma.h"
 
-#undef NAN
-#define NAN octave_NaN
 #define INFINITE lo_ieee_isinf
 #define RUNI oct_randu()
 #define RNOR oct_randn()
@@ -105,7 +103,7 @@
   if (a <= 0 || INFINITE(a))
     {
       for (i=0; i < n; i++)
-        r[i] = NAN;
+        r[i] = octave_numeric_limits<double>::NaN ();
       return;
     }
 
@@ -141,11 +139,9 @@
   return ret;
 }
 
-#undef NAN
 #undef RUNI
 #undef RNOR
 #undef REXP
-#define NAN octave_Float_NaN
 #define RUNI oct_float_randu()
 #define RNOR oct_float_randn()
 #define REXP oct_float_rande()
@@ -162,7 +158,7 @@
   if (a <= 0 || INFINITE(a))
     {
       for (i=0; i < n; i++)
-        r[i] = NAN;
+        r[i] = octave_numeric_limits<float>::NaN ();
       return;
     }
 
--- a/liboctave/numeric/randpoisson.cc	Mon May 16 23:20:50 2016 -0400
+++ b/liboctave/numeric/randpoisson.cc	Tue May 17 01:11:27 2016 -0400
@@ -42,8 +42,6 @@
 #include "randmtzig.h"
 #include "randpoisson.h"
 
-#undef NAN
-#define NAN octave_NaN
 #undef INFINITE
 #define INFINITE lo_ieee_isinf
 #define RUNI oct_randu()
@@ -485,7 +483,7 @@
   if (L < 0.0 || INFINITE(L))
     {
       for (i=0; i<n; i++)
-        p[i] = NAN;
+        p[i] = octave_numeric_limits<double>::NaN ();
     }
   else if (L <= 10.0)
     {
@@ -514,7 +512,7 @@
 oct_randp (double L)
 {
   double ret;
-  if (L < 0.0) ret = NAN;
+  if (L < 0.0) ret = octave_numeric_limits<double>::NaN ();
   else if (L <= 12.0)
     {
       /* From Press, et al. Numerical recipes */
@@ -537,7 +535,7 @@
     {
       /* FIXME: R uses NaN, but the normal approximation suggests that
        * limit should be Inf.  Which is correct? */
-      ret = NAN;
+      ret = octave_numeric_limits<double>::NaN ();
     }
   else
     {
@@ -557,7 +555,7 @@
   if (L < 0.0 || INFINITE(L))
     {
       for (i=0; i<n; i++)
-        p[i] = NAN;
+        p[i] = octave_numeric_limits<double>::NaN ();
     }
   else if (L <= 10.0)
     {
@@ -587,7 +585,7 @@
 {
   double L = FL;
   float ret;
-  if (L < 0.0) ret = NAN;
+  if (L < 0.0) ret = octave_numeric_limits<float>::NaN ();
   else if (L <= 12.0)
     {
       /* From Press, et al. Numerical recipes */
@@ -610,7 +608,7 @@
     {
       /* FIXME: R uses NaN, but the normal approximation suggests that
        * limit should be Inf. Which is correct? */
-      ret = NAN;
+      ret = octave_numeric_limits<float>::NaN ();
     }
   else
     {
--- a/liboctave/util/lo-ieee.h	Mon May 16 23:20:50 2016 -0400
+++ b/liboctave/util/lo-ieee.h	Tue May 17 01:11:27 2016 -0400
@@ -121,18 +121,24 @@
 struct octave_numeric_limits
 {
   static T NA (void) { return static_cast<T> (0); }
+  static T NaN (void) { return static_cast<T> (0); }
+  static T Inf (void) { return static_cast<T> (0); }
 };
 
 template <>
 struct octave_numeric_limits<double>
 {
   static double NA (void) { return octave_NA; }
+  static double NaN (void) { return octave_NaN; }
+  static double Inf (void) { return octave_Inf; }
 };
 
 template <>
 struct octave_numeric_limits<float>
 {
   static float NA (void) { return octave_Float_NA; }
+  static float NaN (void) { return octave_Float_NaN; }
+  static float Inf (void) { return octave_Float_Inf; }
 };
 
 #endif