optimize(); replace -Inf by *negative*; warning msg w/ more detail wh…

…en replacing NA or +/-Inf git-svn-id: https://svn.r-project.org/R/trunk@87693 00db46b3-68df-0310-9c12-caf00c1e9a41
r-devel · Feb 6, 2025 · c24cbaa · c24cbaa
1 parent 82fc32e
commit c24cbaa
Show file tree

Hide file tree

Showing 3 changed files with 56 additions and 16 deletions.
diff --git a/doc/NEWS.Rd b/doc/NEWS.Rd
@@ -191,6 +191,10 @@
 
       \item \code{unique()}'s default method now also deals with
       \code{"difftime"} objects.
+
+      \item \code{optimize(f, *)} when \code{f(x)} is not finite says
+      more about the value in its \code{warning} message.  It no longer
+      replaces \code{-Inf} by the largest \emph{positive} finite number.
     }
   }
 

diff --git a/src/library/stats/src/optimize.c b/src/library/stats/src/optimize.c
@@ -1,8 +1,8 @@
 /*
  *  R : A Computer Language for Statistical Data Analysis
+ *  Copyright (C) 1998--2025  The R Core Team
  *  Copyright (C) 1995, 1996  Robert Gentleman and Ross Ihaka
  *  Copyright (C) 2003-2004  The R Foundation
- *  Copyright (C) 1998--2023  The R Core Team
  *
  *  This program is free software; you can redistribute it and/or modify
  *  it under the terms of the GNU General Public License as published by
@@ -219,8 +219,13 @@ static double fcn1(double x, void *arg_info)
     case REALSXP:
 	if (length(s) != 1) goto badvalue;
 	if (!R_FINITE(REAL(s)[0])) {
-	    warning(_("NA/Inf replaced by maximum positive value"));
-	    return DBL_MAX;
+	    if(REAL(s)[0] == R_NegInf) { // keep sign for root finding !
+		warning(_("-Inf replaced by maximally negative value"));
+		return -DBL_MAX;
+	    } else {
+		warning(_("%s replaced by maximum positive value"), ISNAN(REAL(s)[0]) ? "NA/NaN" : "Inf");
+		return DBL_MAX;
+	    }
 	}
 	else return REAL(s)[0];
 	break;
@@ -232,33 +237,31 @@ static double fcn1(double x, void *arg_info)
     return 0;/* for -Wall */
 }
 
-/* fmin(f, xmin, xmax tol) */
+/* Called from optimize() as
+ * .External2(C_do_fmin,  function(arg) +/- f(arg, ...), lower, upper, tol)
+ * fmin(f, xmin, xmax tol) */
 SEXP do_fmin(SEXP call, SEXP op, SEXP args, SEXP rho)
 {
-    double xmin, xmax, tol;
-    SEXP v, res;
-    struct callinfo info;
-
     args = CDR(args);
     PrintDefaults();
 
     /* the function to be minimized */
 
-    v = CAR(args);
+    SEXP v = CAR(args);
     if (!isFunction(v))
 	error(_("attempt to minimize non-function"));
     args = CDR(args);
 
     /* xmin */
 
-    xmin = asReal(CAR(args));
+    double xmin = asReal(CAR(args));
     if (!R_FINITE(xmin))
 	error(_("invalid '%s' value"), "xmin");
     args = CDR(args);
 
     /* xmax */
 
-    xmax = asReal(CAR(args));
+    double xmax = asReal(CAR(args));
     if (!R_FINITE(xmax))
 	error(_("invalid '%s' value"), "xmax");
     if (xmin >= xmax)
@@ -267,13 +270,14 @@ SEXP do_fmin(SEXP call, SEXP op, SEXP args, SEXP rho)
 
     /* tol */
 
-    tol = asReal(CAR(args));
+    double tol = asReal(CAR(args));
     if (!R_FINITE(tol) || tol <= 0.0)
 	error(_("invalid '%s' value"), "tol");
 
+    struct callinfo info;
     info.R_env = rho;
     PROTECT(info.R_fcall = lang2(v, R_NilValue));
-    PROTECT(res = allocVector(REALSXP, 1));
+    SEXP res = PROTECT(allocVector(REALSXP, 1));
     REAL(res)[0] = Brent_fmin(xmin, xmax, fcn1, &info, tol);
     UNPROTECT(2);
     return res;
@@ -309,7 +313,7 @@ static double fcn2(double x, void *arg_info)
 		warning(_("-Inf replaced by maximally negative value"));
 		return -DBL_MAX;
 	    } else {
-		warning(_("NA/Inf replaced by maximum positive value"));
+		warning(_("%s replaced by maximum positive value"), ISNAN(REAL(s)[0]) ? "NA/NaN" : "Inf");
 		return DBL_MAX;
 	    }
 	}
@@ -527,8 +531,13 @@ static void fcn(int n, double *x, double *f, void *arg_state)
     case REALSXP:
 	if (length(s) != 1) goto badvalue;
 	if (!R_FINITE(REAL(s)[0])) {
-	    warning(_("NA/Inf replaced by maximum positive value"));
-	    *f = DBL_MAX;
+	    if(REAL(s)[0] == R_NegInf) { // keep sign for root finding !
+		warning(_("-Inf replaced by maximally negative value"));
+		*f = -DBL_MAX;
+	    } else {
+		warning(_("%s replaced by maximum positive value"), ISNAN(REAL(s)[0]) ? "NA/NaN" : "Inf");
+		*f = DBL_MAX;
+	    }
 	}
 	else *f = REAL(s)[0];
 	break;

diff --git a/tests/reg-tests-1e.R b/tests/reg-tests-1e.R
@@ -1806,6 +1806,33 @@ stopifnot(inherits(unidt, "difftime"), length(unidt) <= 2) # '2': allow "inaccur
 ## unique() lost the class  in R < 4.5.0
 
 
+## optimize(f(x), *) message when f(x) is not finite
+ff <- function(x) ifelse(x < -10, (x+10)*exp(x^2),
+                     ifelse(x > 100, NaN,
+                        ifelse(x > 30, exp((x-20)^2),
+                               (4 - x)^2)))
+cf <- as.data.frame(curve(ff, -20, 120, ylim = c(-2,200)))
+str(ok <- optimize(ff, c(-10, 10)))
+stopifnot(all.equal(list(minimum = 4, objective = 0), ok))
+op <- options(warn=0)
+str(of2 <- optimize(ff, c(-140, 250))); summary(warnings()); uw2 <- unique(warnings())
+## NA/NaN  and -Inf  (no +Inf)
+str(of3 <- optimize(ff, c(-20, 120)));  summary(warnings()); uw3 <- unique(warnings())
+## only 1 Inf
+str(of4 <- optimize(ff, c(-10, 180)));  summary(warnings()); uw4 <- unique(warnings())
+## +Inf and many NA/NaN
+c(uw2, uw3, uw4)
+stopifnot(all.equal(of3, ok),
+          identical(c(2:1,2L), lengths(list(uw2, uw3, uw4))))
+if(englishMsgs)
+    stopifnot(identical(c("-Inf replaced by maximally negative value",
+                           "Inf replaced by maximum positive value",
+                        "NA/NaN replaced by maximum positive value"),
+                        sort(unique(c(names(uw2), names(uw3), names(uw4))))))
+options(op)# reverting
+## in R < 4.4.z  only *one* message .. "NA/Inf replaced by ...."
+
+
 
 ## keep at end
 rbind(last =  proc.time() - .pt,