c(NA, 1i) also coerces to cpl(Re=NA, Im=0) ; ditto cumsum() etc

git-svn-id: https://svn.r-project.org/R/trunk@85525 00db46b3-68df-0310-9c12-caf00c1e9a41

doc/NEWS.Rd

@@ -212,6 +212,11 @@
 
       \item \code{as.complex("1i")} now returns \code{1i} instead of
       \code{NA} with a warning.
+
+      \item \code{z <- c(NA, 1i)} now keeps the imaginary part
+      \code{Im(z[1]) == 0}, no longer coercing to \code{NA_complex_}.
+      Simlarly, \code{cumsum(z)} correctly sums real and imaginary parts
+      separately, i.e., without \dQuote{crosstalk} in case of \code{NA}s.
     }
   }
 
@@ -338,7 +343,7 @@
       \code{xfig()} are now deprecated.  They do not support recent
       graphics enhancements and their font-handling is rudimentary. The
       intention is to retain them for historical interest as long as
-      they remain somewhat functional.   
+      they remain somewhat functional.
     }
   }
 
@@ -465,7 +470,7 @@
       \PR{18612} (and more) reported by Mikael Jagan.
 
       \item \code{list.files()} on Windows now returns also files with names
-      longer that 260 bytes (the Windows limit is 260 characters). 
+      longer that 260 bytes (the Windows limit is 260 characters).
       Previously, some file names particularly with \sQuote{East Asian}
       characters were omitted.
 
@@ -484,12 +489,13 @@
       \PR{17780}).
 
       \item Datetime functions now restore even an empty \env{TZ}
-      environment variable after internal time zone changes (\PR{17724}). 
+      environment variable after internal time zone changes (\PR{17724}).
       This makes results of datetime functions with this (typically
       unintentional) setting more predictable.
     }
   }
 }
+
 \section{\Rlogo CHANGES IN R 4.3.2 patched}{
   \subsection{BUG FIXES}{
     \itemize{

src/library/base/man/c.Rd

@@ -1,6 +1,6 @@
 % File src/library/base/man/c.Rd
 % Part of the R package, https://www.R-project.org
-% Copyright 1995-2022 R Core Team
+% Copyright 1995-2023 R Core Team
 % Distributed under GPL 2 or later
 \name{c}
 \title{Combine Values into a Vector or List}
@@ -37,7 +37,11 @@ c(\dots)
   are treated as one-element \code{\link{list}}s
   which cannot be unlisted even if \code{recursive = TRUE}.
 
-    There is a \code{\link{c.factor}} method which combines factors into
+  If the output type is \code{\link{complex}}, logical, integer, and double
+  \code{NA}s keep their imaginary parts zero when coerced, and hence will
+  \emph{not} become \code{NA_complex_} (with imaginary part \code{NA}).
+
+  There is a \code{\link{c.factor}} method which combines factors into
   a factor.
 
   \code{c} is sometimes used for its side effect of removing attributes

src/library/base/man/cumsum.Rd

@@ -1,6 +1,6 @@
 % File src/library/base/man/cumsum.Rd
 % Part of the R package, https://www.R-project.org
-% Copyright 1995-2009 R Core Team
+% Copyright 1995-2023 R Core Team
 % Distributed under GPL 2 or later
 
 \name{cumsum}
@@ -35,6 +35,9 @@ cummin(x)
   An \code{NA} value in \code{x} causes the corresponding and following
   elements of the return value to be \code{NA}, as does integer overflow
   in \code{cumsum} (with a warning).
+  In the complex case with \code{\link{NA}}s, these \code{NA} elements may
+  have finite real or imaginary parts, notably for \code{cumsum()},
+  fulfulling the identity \code{Im(cumsum(x))} \eqn{\equiv}{=id=} \code{cumsum(Im(x))}.
 }
 \section{S4 methods}{
   \code{cumsum} and \code{cumprod} are S4 generic functions:

src/main/bind.c

@@ -417,7 +417,11 @@ ComplexAnswer(SEXP x, struct BindData *data, SEXP call)
 	    xi = LOGICAL(x)[i];
 	    if (xi == NA_LOGICAL) {
 		COMPLEX(data->ans_ptr)[data->ans_length].r = NA_REAL;
+#ifdef NA_TO_COMPLEX_NA
 		COMPLEX(data->ans_ptr)[data->ans_length].i = NA_REAL;
+#else
+		COMPLEX(data->ans_ptr)[data->ans_length].i = 0.0;
+#endif
 	    }
 	    else {
 		COMPLEX(data->ans_ptr)[data->ans_length].r = xi;
@@ -431,7 +435,11 @@ ComplexAnswer(SEXP x, struct BindData *data, SEXP call)
 	    xi = INTEGER(x)[i];
 	    if (xi == NA_INTEGER) {
 		COMPLEX(data->ans_ptr)[data->ans_length].r = NA_REAL;
+#ifdef NA_TO_COMPLEX_NA
 		COMPLEX(data->ans_ptr)[data->ans_length].i = NA_REAL;
+#else
+		COMPLEX(data->ans_ptr)[data->ans_length].i = 0.0;
+#endif
 	    }
 	    else {
 		COMPLEX(data->ans_ptr)[data->ans_length].r = xi;

src/main/cum.c

@@ -1,7 +1,7 @@
 /*
  *  R : A Computer Language for Statistical Data Analysis
+ *  Copyright (C) 1997--2023  The R Core Team
  *  Copyright (C) 1995, 1996  Robert Gentleman and Ross Ihaka
- *  Copyright (C) 1997--2022  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
@@ -73,22 +73,23 @@ static SEXP icumsum(SEXP x, SEXP s)
     return s;
 }
 
-static SEXP chandleNaN(SEXP x, SEXP s)
+/* For complex result: recompute once we know one of the result's {re, im} fulfills  ISNAN(.),
+   (speed optimized for the case of *no* NA|NaN) : */
+static SEXP chandleNaN(SEXP x, SEXP s, Rboolean r_isN, Rboolean i_isN)
 {
     Rboolean hasNA = FALSE;
     Rboolean hasNaN = FALSE;
 
     for (R_xlen_t i = 0 ; i < XLENGTH(x) ; i++) {
 	hasNaN = hasNaN || ISNAN(COMPLEX(x)[i].r) || ISNAN(COMPLEX(x)[i].i);
-	hasNA = hasNA || (hasNaN && (R_IsNA(COMPLEX(x)[i].r)
-	                             || R_IsNA(COMPLEX(x)[i].i)));
-
+	hasNA = hasNA || (hasNaN && (R_IsNA(COMPLEX(x)[i].r) ||
+				     R_IsNA(COMPLEX(x)[i].i)));
 	if (hasNA) {
-	    COMPLEX(s)[i].r = NA_REAL;
-	    COMPLEX(s)[i].i = NA_REAL;
+	    if(r_isN) COMPLEX(s)[i].r = NA_REAL;
+	    if(i_isN) COMPLEX(s)[i].i = NA_REAL;
 	} else if (hasNaN) {
-	    COMPLEX(s)[i].r = R_NaN;
-	    COMPLEX(s)[i].i = R_NaN;
+	    if(r_isN) COMPLEX(s)[i].r = R_NaN;
+	    if(i_isN) COMPLEX(s)[i].i = R_NaN;
 	}
     }
     return s;
@@ -105,7 +106,7 @@ static SEXP ccumsum(SEXP x, SEXP s)
 	COMPLEX(s)[i].r = sum.r;
 	COMPLEX(s)[i].i = sum.i;
     }
-    return (ISNAN(sum.r) || ISNAN(sum.i)) ? chandleNaN(x, s) : s;
+    return (ISNAN(sum.r) || ISNAN(sum.i)) ? chandleNaN(x, s, ISNAN(sum.r), ISNAN(sum.i)) : s;
 }
 
 static SEXP cumprod(SEXP x, SEXP s)
@@ -133,7 +134,7 @@ static SEXP ccumprod(SEXP x, SEXP s)
 	COMPLEX(s)[i].r = prod.r;
 	COMPLEX(s)[i].i = prod.i;
     }
-    return (ISNAN(prod.r) || ISNAN(prod.i)) ? chandleNaN(x, s) : s;
+    return (ISNAN(prod.r) || ISNAN(prod.i)) ? chandleNaN(x, s, ISNAN(prod.r), ISNAN(prod.i)) : s;
 }
 
 static SEXP cummax(SEXP x, SEXP s)

tests/reg-tests-1a.R

@@ -3066,10 +3066,24 @@ stopifnot(identical(cumprod(x), r))
 stopifnot(identical(cummin(x), r))
 stopifnot(identical(cummax(x), r))
 # complex
-x <- c(1+1i, NA, 3)
-r <- c(1+1i, NA, NA)
-stopifnot(identical(cumsum(x), r))
-stopifnot(identical(cumprod(x), r))
+cx <- function(r,i) complex(real=r, imaginary=i)
+NA.1 <- cx(NA, 1)
+NA.2 <- cx(NA, 2)
+NA_C <- NA_complex_ # = complex(r=NA, i=NA)
+y <- x <- c(1+1i, NA, 3)
+stopifnot(identical(x[2], cx(NA,0))) # newly true
+y[2] <- NA.1
+stopifnot(exprs = {
+    identical(Im(cumsum(x)), cumsum(Im(x)))
+    identical(Re(cumsum(x)), cumsum(Re(x)))
+    identical(Im(cumsum(y)), cumsum(Im(y)))
+    identical(Re(cumsum(y)), cumsum(Re(y)))
+    identical( sum(x), tail( cumsum(x), 1L))
+    identical(prod(x), tail(cumprod(x), 1L))
+    identical(cumsum (x), c(1+1i, NA.1, NA.1)) # new
+    identical(cumsum (y), c(1+1i, NA.2, NA.2)) #  "
+    identical(cumprod(x), c(1+1i, NA_C, NA_C))
+})
 # integer
 x <- c(1L, NA, 3L)
 r <- c(1L, NA, NA)

tests/reg-tests-1e.R

@@ -1048,6 +1048,25 @@ stopifnot(exprs = {
 ## returned NA_complex_ *with* a warning, in R <= 4.4.0
 
 
+## c(NA, <cplx>); cumsum(<cplx_w_NA>) -- related to as.complex(NA_real_); R-devel ML
+cx <- function(r,i) complex(real=r, imaginary=i)
+pz <- function(z) noquote(paste0("(", Re(z), ",", Im(z), ")"))
+pz(z <- c(1i, NA))
+NA.1 <- cx(NA, 1)
+stopifnot(exprs = {
+    Im(z) == 1:0 # was  (TRUE  NA)
+    identical(z, cx(c(0,NA), c(1,0))) # new
+    identical(sum(z), NA.1)           #  "
+    ## these were all TRUE already :
+    identical(prod(z), NA_complex_)
+    identical(Im(cumsum(z)), cumsum(Im(z)))
+    identical(Re(cumsum(z)), cumsum(Re(z)))
+    identical( sum(z), tail( cumsum(z), 1L))
+    identical(prod(z), tail(cumprod(z), 1L))
+})
+##
+
+
 
 ## keep at end
 rbind(last =  proc.time() - .pt,