Bug fixes in trivial cases of sumsymouterSparse

baddstats · baddstats · commit 7c601ff3e721 · 2026-04-24T14:06:11.000+08:00
diff --git a/DESCRIPTION b/DESCRIPTION
@@ -1,6 +1,6 @@
 Package: spatstat.sparse
-Version: 3.1-0.007
-Date: 2026-04-23
+Version: 3.1-0.009
+Date: 2026-04-24
 Title: Sparse Three-Dimensional Arrays and Linear Algebra Utilities
 Authors@R: c(person("Adrian", "Baddeley", 
                     role = c("aut", "cre", "cph"),
diff --git a/NAMESPACE b/NAMESPACE
@@ -42,6 +42,7 @@ export("evalSparse3Dentrywise")
 export("expandSparse")
 export("gridadjacencymatrix")
 export("inside3Darray")
+export("length.sparse3Darray")
 export("mapSparseEntries")
 export("marginSumsSparse")
 export("Math.sparse3Darray")
@@ -84,6 +85,7 @@ S3method("Summary", "sparse3Darray")
 S3method("as.array", "sparse3Darray")
  S3method("dimnames", "sparse3Darray")
  S3method("dim", "sparse3Darray")
+ S3method("length", "sparse3Darray")
  S3method("print", "sparse3Darray")
   S3method("[", "sparse3Darray")
 # ......................................... 
diff --git a/NEWS b/NEWS
@@ -1,12 +1,16 @@
 
-        CHANGES IN spatstat.sparse VERSION 3.1-0.007
+        CHANGES IN spatstat.sparse VERSION 3.1-0.009
 
 OVERVIEW
 
     o Simulate a Markov chain with sparse transition matrix.
 
     o Recreate a symmetric matrix from its lower triangular entries.
 
+    o Weighted sum of outer products in a 3D array or sparse array.
+
+    o Minor improvements.
+
 NEW FUNCTIONS
 
     o runSparseMarkovChain
@@ -17,6 +21,12 @@ NEW FUNCTIONS
     Recreate a symmetric matrix from a vector containing its lower triangular
     (or upper triangular) entries.
 
+    o length.sparse3Darray
+    Method for 'length' for sparse 3D arrays.
+
+    o sumsymouter
+    Weighted sum of outer products in a 3D array or sparse array
+
         CHANGES IN spatstat.sparse VERSION 3.1-0
 
 OVERVIEW
diff --git a/R/linalg.R b/R/linalg.R
@@ -6,7 +6,7 @@
 #'  Copyright (c) Adrian Baddeley, Ege Rubak and Rolf Turner 2016-2020
 #'  GNU Public Licence >= 2.0
 #'
-#' $Revision: 1.38 $ $Date: 2022/05/23 02:33:06 $
+#' $Revision: 1.39 $ $Date: 2026/04/24 05:57:32 $
 #'
 
 sumouter <- function(x, w=NULL, y=x) {
@@ -239,8 +239,8 @@ bilinearform <- function(x, v, y) {
 }
 
 sumsymouter <- function(x, w=NULL, distinct=TRUE) {
-  ## x is a 3D array
-  ## w is a matrix
+  ## x is a 3D array or sparse 3D array
+  ## w is a matrix or sparse matrix
   ## Computes the sum of outer(x[,i,j], x[,j,i]) * w[i,j] over all pairs i != j
   ## handle complex values
   if(is.complex(w)) {
@@ -256,20 +256,25 @@ sumsymouter <- function(x, w=NULL, distinct=TRUE) {
     result <- a - b + (d - a - b) * 1i
     return(result)
   }
-  ## handle sparse arrays
-  if(inherits(x, c("sparseSlab", "sparse3Darray")) &&
-     (is.null(w) || inherits(w, "sparseMatrix")))
-    return(sumsymouterSparse(x, w, distinct=distinct))
   ## arguments are numeric
-  x <- as.array(x)
+  ## validate dimensions
   stopifnot(length(dim(x)) == 3)
   if(dim(x)[2L] != dim(x)[3L])
     stop("The second and third dimensions of x should be equal")
   if(!is.null(w)) {
-    w <- as.matrix(w)
+    stopifnot(length(dim(w)) == 2)
     if(!all(dim(w) == dim(x)[-1L]))
       stop("Dimensions of w should match the second and third dimensions of x")
   }
+  ## handle sparse arrays
+  if(inherits(x, c("sparseSlab", "sparse3Darray"))) {
+    if(!is.null(w))
+      w <- as(w, "sparseMatrix")
+    return(sumsymouterSparse(x, w, distinct=distinct))
+  }
+  x <- as.array(x)
+  if(!is.null(w)) 
+    w <- as.matrix(w)
   p <- dim(x)[1L]
   n <- dim(x)[2L]
   if(!distinct) {
diff --git a/R/sparse3Darray.R b/R/sparse3Darray.R
@@ -126,6 +126,8 @@ as.sparse3Darray <- function(x, ...) {
   return(y)
 }
 
+length.sparse3Darray <- function(x) { prod(dim(x)) }
+
 dim.sparse3Darray <- function(x) { x$dim }
 
 "dim<-.sparse3Darray" <- function(x, value) {
diff --git a/R/sparselinalg.R b/R/sparselinalg.R
@@ -6,7 +6,7 @@
 #'   Copyright (c) Adrian Baddeley, Ege Rubak and Rolf Turner 2016-2020
 #'   GNU Public Licence >= 2.0
 #' 
-#'   $Revision: 1.19 $  $Date: 2021/01/08 01:16:41 $
+#'   $Revision: 1.21 $  $Date: 2026/04/24 05:22:07 $
 
 marginSumsSparse <- function(X, MARGIN) {
   #' equivalent to apply(X, MARGIN, sum)
@@ -221,17 +221,17 @@ sumsymouterSparse <- function(x, w=NULL, distinct=TRUE, dbg=FALSE) {
   if(distinct) {
     #' remove entries with j = k
     ok <- with(df, j != k)
-    df <- df[ok, , drop=TRUE]
+    df <- df[ok, , drop=FALSE]
   }
   ## trivial?
-  if(nrow(df) < 2) {
+  if(nrow(df) == 0) {
     y <- matrix(0, m, m)
     dimnames(y) <- rep(dimnames(x)[1], 2)
     return(y)
   }
   ## order by increasing j, then k
   oo <- with(df, order(j, k, i))
-  df <- df[oo, ]
+  df <- df[oo, , drop=FALSE]
   ## now provide ordering by increasing k then j
   ff <- with(df, order(k,j,i))
   ##
diff --git a/inst/doc/packagesizes.txt b/inst/doc/packagesizes.txt
@@ -19,4 +19,4 @@ date version nhelpfiles nobjects ndatasets Rlines srclines
 "2023-03-12" "3.0-1" 15 48 0 2092 740
 "2023-10-24" "3.0-3" 15 48 0 2092 740
 "2024-06-21" "3.1-0" 15 48 0 2092 740
-"2026-04-23" "3.1-0.007" 17 50 0 2217 943
+"2026-04-24" "3.1-0.009" 18 51 0 2224 943
diff --git a/inst/info/packagesizes.txt b/inst/info/packagesizes.txt
@@ -19,4 +19,4 @@ date version nhelpfiles nobjects ndatasets Rlines srclines
 "2023-03-12" "3.0-1" 15 48 0 2092 740
 "2023-10-24" "3.0-3" 15 48 0 2092 740
 "2024-06-21" "3.1-0" 15 48 0 2092 740
-"2026-04-23" "3.1-0.007" 17 50 0 2217 943
+"2026-04-24" "3.1-0.009" 18 51 0 2224 943
diff --git a/man/methods.sparse3Darray.Rd b/man/methods.sparse3Darray.Rd
@@ -6,6 +6,7 @@
 \alias{dim<-.sparse3Darray}
 \alias{dimnames.sparse3Darray}
 \alias{dimnames<-.sparse3Darray}
+\alias{length.sparse3Darray}
 \alias{print.sparse3Darray}
 \title{
   Methods for Sparse Three-Dimensional Arrays
@@ -20,6 +21,7 @@
  \method{dim}{sparse3Darray}(x) <- value
  \method{dimnames}{sparse3Darray}(x)
  \method{dimnames}{sparse3Darray}(x) <- value
+ \method{length}{sparse3Darray}(x)
  \method{print}{sparse3Darray}(x, \dots)
 }
 \arguments{
@@ -40,8 +42,9 @@
   \code{\link[base]{dim<-}},
   \code{\link[base]{dimnames}},
   \code{\link[base]{dimnames<-}}
+  \code{\link[base]{length}}
   and \code{\link[base]{print}}
-  for the class \code{"sparse#Darray"} of sparse three-dimensional
+  for the class \code{"sparse3Darray"} of sparse three-dimensional
   arrays.
 
   For \code{dimnames(x) <- value}, the \code{value} should either be
@@ -52,7 +55,10 @@
   vector of length 3 giving the new dimensions of the array. Note that this
   operation does not change the array positions of the non-zero entries
   (unlike \code{dim(x) <- value} for a full array). An error occurs if
-  some of the non-zero entries would lie outside the new extent of the array.
+  some of the non-zero entries would lie outside the new extent of the
+  array.
+
+  The \code{length} of a sparse array is the product of its dimensions.
 }
 \value{
   \code{anyNA} returns a single logical value.
@@ -64,6 +70,8 @@
 
   \code{dim<-} and \code{dimnames<-} return a sparse 3D array.
 
+  \code{length} returns an integer.
+  
   \code{print} returns \code{NULL}, invisibly.
 }
 \author{
diff --git a/man/spatstat.sparse-internal.Rd b/man/spatstat.sparse-internal.Rd
@@ -6,7 +6,6 @@
 %%%%  undocumented linear algebra
 \alias{checksolve}
 \alias{check.mat.mul}
-\alias{sumsymouter}
 %% sparse 3D arrays
 \alias{unionOfSparseIndices}
 \alias{inside3Darray}
@@ -31,7 +30,6 @@ representativeRows(x)
 %%% undocumented linear algebra
 checksolve(M, action, descrip, target)
 check.mat.mul(A, B, Acols, Brows, fatal)
-sumsymouter(x, w, distinct)
 %% sparse 3D arrays
 unionOfSparseIndices(A,B)
 inside3Darray(d, i, j, k)
diff --git a/man/sumouter.Rd b/man/sumouter.Rd
@@ -105,5 +105,6 @@
   and \rolf
 }
 \keyword{array}
+\keyword{algebra}
  
  
diff --git a/man/sumsymouter.Rd b/man/sumsymouter.Rd
@@ -0,0 +1,70 @@
+\name{sumsymouter}
+\alias{sumsymouter}
+\title{Compute Quadratic Forms in 3D Arrays}
+\description{
+  Calculates certain quadratic forms of 3D arrays.
+}
+\usage{
+  sumsymouter(x, w = NULL, distinct = TRUE)
+}
+\arguments{
+  \item{x}{
+    Three-dimensional array (class \code{array})
+    or sparse three-dimensional array (class \code{sparse3Darray}).
+  }
+  \item{w}{Matrix or sparse matrix.}
+  \item{distinct}{
+    Logical value indicating whether to exclude diagonal contributions.
+    See Details.
+  }
+}
+\details{
+  Given a three-dimensional matrix or sparse matrix \eqn{x}
+  and a matrix or sparse matrix \eqn{w},
+  this function computes the matrix \eqn{M} with entries
+  \deqn{
+    M_{i,j} = \sum_k \sum_m w_{k,m} x_{i, k, m} x_{j, m, k}
+  }{
+    M[i,j] = sum[k] sum[m] w[k,m] * x[i,k,m] * x[j, m, k]
+  },
+  the weighted sum of symmetric outer products of \eqn{x}
+  with weights \eqn{w}.
+
+  If \code{distinct=TRUE} (the default), terms with \eqn{k=m} are
+  omitted from the sum. If \code{distinct=FALSE} they are included.
+
+  The array \code{x} should have dimensions
+  \eqn{p \times n \times n}{p * n * n}
+  and the matrix \code{w} should have dimensions \eqn{n \times n}{n * n}
+  if it is given. 
+  The result \eqn{M} is a \eqn{p \times p}{p * p} square matrix.
+
+  If \code{w} is missing or \code{NULL}, all weights are taken to be 1.
+
+  Entries in \code{x} and \code{w} may be integer, numeric, complex
+  or logical values.
+  
+  The result is a symmetric matrix if \code{w} is a symmetric matrix,
+  or if \code{w} is missing or \code{NULL}. Otherwise it is not
+  guaranteed to be symmetric.
+}
+\value{
+  A square matrix with numeric or complex values.
+}
+\author{
+  \adrian.
+}
+\seealso{
+   \code{\link{sumouter}}
+}
+\examples{
+   d <- c(2, 3, 3)
+   n <- prod(d)
+   v <- ifelse(runif(n) < 0.6, 0, round(runif(n), 1))
+   x <- array(v, dim=d)
+   w <- outer(1:3, 1:3)
+   sumsymouter(x, w)
+}
+\keyword{array}
+\keyword{algebra}
+\concept{sparse}
diff --git a/src/spasumsymout.h b/src/spasumsymout.h
@@ -13,7 +13,7 @@
 
   WEIGHTS   (#ifdef) use weights 
 
-  $Revision: 1.11 $  $Date: 2022/10/19 08:16:51 $
+  $Revision: 1.12 $  $Date: 2026/04/24 05:36:26 $
 
   Copyright (C) Adrian Baddeley, Ege Rubak and Rolf Turner 2001-2018
   Licence: GNU Public Licence >= 2
@@ -62,7 +62,7 @@ void FNAME(
   R = *lenw;
 #endif
 
-  if(L <= 1 || N <= 1 || M <= 0) return;
+  if(L <= 0 || N <= 1 || M <= 0) return;
 
   /* Create space to store array in k-major order*/
   it = (int *) R_alloc(L, sizeof(int));

Original file line number	Diff line number	Diff line change
`@@ -126,6 +126,8 @@ as.sparse3Darray <- function(x, ...) {`
`126`	`126`	`return(y)`
`127`	`127`	`}`
`128`	`128`
	`129`	`+length.sparse3Darray <- function(x) { prod(dim(x)) }`
	`130`	`+`
`129`	`131`	`dim.sparse3Darray <- function(x) { x$dim }`
`130`	`132`
`131`	`133`	`"dim<-.sparse3Darray" <- function(x, value) {`