R Under development (unstable) (2024-03-19 r86153 ucrt) -- "Unsuffered Consequences"
Copyright (C) 2024 The R Foundation for Statistical Computing
Platform: x86_64-w64-mingw32/x64

R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.

R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.

Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.

> library(Matrix)
> ### Do all kinds of object creation and coercion
> 
> source(system.file("test-tools.R", package = "Matrix"))
Loading required package: tools
> 
> .drop.factors <- function(x, check = FALSE)
+     `slot<-`(x, "factors", check = check, value = list())
> 
> ## the empty ones:
> checkMatrix(new("dgeMatrix"))
 norm(m [0 x 0]) : 1  I  F  M  ok
 Summary:  ok
2*m =?= m+m: identical
m >= m for all: ok
symmpart(m) + skewpart(m) == m: ok;  determinant(): ok
> checkMatrix(Matrix(,0,0))
 norm(m [0 x 0]) : 1  I  F  M  ok
 Summary:  ok
2*m =?= m+m: identical
m >= m for all: ok
symmpart(m) + skewpart(m) == m: ok;  determinant(): ok
> 
> ## "dge"
> assertError( new("dgeMatrix", Dim = c(2,2), x= 1:4) )# double 'Dim'
> assertError( new("dgeMatrix", Dim = as.integer(c(2,2)), x= 1:4) )# int 'x'
> assertError( new("dgeMatrix", Dim = 2:2, x=as.double(1:4)) )# length(Dim) !=2
> assertError( new("dgeMatrix", Dim = as.integer(c(2,2)), x= as.double(1:5)))
> 
> checkMatrix(m1 <- Matrix(1:6, ncol=2))
 norm(m [3 x 2]) : 1  I  F  M  ok
 Summary:  ok
2*m =?= m+m: identical
m >= m for all: ok
m < m for none: ok
> checkMatrix(m2 <- Matrix(1:7 +0, ncol=3)) # a (desired) warning
 norm(m [3 x 3]) : 1  I  F  M  ok
 Summary:  ok
2*m =?= m+m: identical
m >= m for all: ok
m < m for none: ok
symmpart(m) + skewpart(m) == m: ok;  determinant(): ok
Warning message:
In Matrix(1:7 + 0, ncol = 3) :
  data length [7] is not a sub-multiple or multiple of the number of rows [3]
> c("dgeMatrix", "ddenseMatrix", "generalMatrix", "dMatrix",
+   "denseMatrix", "Matrix") -> m1.cl
> stopifnot(!anyNA(match(m1.cl, is(m1))),
+ 	  dim(t(m1)) == 2:3, identical(m1, t(t(m1))))
> c.nam <- paste("C",1:2, sep='')
> dimnames(m1) <- list(NULL, c.nam)
> checkMatrix(m1) # failed in 0.999375-10
 norm(m [3 x 2]) : 1  I  F  M  ok
 Summary:  ok
2*m =?= m+m: identical
m >= m for all: ok
m < m for none: ok
> checkMatrix(tm1 <- t(m1))
 norm(m [2 x 3]) : 1  I  F  M  ok
 Summary:  ok
2*m =?= m+m: identical
m >= m for all: ok
m < m for none: ok
> stopifnot(colnames(m1) == c.nam,
+ 	  identical(dimnames(tm1), list(c.nam, NULL)),
+ 	  identical(m1, t(tm1)))
> 
> ## an example of *named* dimnames
> (t34N <- as(unclass(table(x = gl(3,4), y=gl(4,3))), "unpackedMatrix"))
3 x 4 Matrix of class "dgeMatrix"
   y
x   1 2 3 4
  1 3 1 0 0
  2 0 2 2 0
  3 0 0 1 3
> stopifnot(identical(dimnames(t34N),
+ 		    dimnames(as(t34N, "matrix"))),
+           identical(t34N, t(t(t34N))))
> 
> ## "dpo"
> checkMatrix(cm <- crossprod(m1))
 norm(m [2 x 2]) : 1  I  F  M  ok
 Summary:  ok
2*m =?= m+m: ok
m >= m for all: ok
m < m for none: ok
symmpart(m) + skewpart(m) == m: ok;  determinant(): ok
> checkMatrix(cp <- as(cm, "packedMatrix"))# 'dpp' + factors
 norm(m [2 x 2]) : 1  I  F  M  ok
 Summary:  ok
2*m =?= m+m: ok
m >= m for all: ok
m < m for none: ok
symmpart(m) + skewpart(m) == m: ok;  determinant(): ok
> checkMatrix(cs <- as(cm, "dsyMatrix"))# 'dsy' + factors
 norm(m [2 x 2]) : 1  I  F  M  ok
 Summary:  ok
2*m =?= m+m: identical
m >= m for all: ok
m < m for none: ok
symmpart(m) + skewpart(m) == m: ok;  determinant(): ok
> checkMatrix(dcm <- as(cm, "generalMatrix"))#'dge'
 norm(m [2 x 2]) : 1  I  F  M  ok
 Summary:  ok
2*m =?= m+m: identical
m >= m for all: ok
m < m for none: ok
symmpart(m) + skewpart(m) == m: ok;  determinant(): ok
> checkMatrix(mcm <- as(cm, "dMatrix")) # 'dsy' + factors -- buglet? rather == cm?
 norm(m [2 x 2]) : 1  I  F  M  ok
 Summary:  ok
2*m =?= m+m: identical
m >= m for all: ok
m < m for none: ok
symmpart(m) + skewpart(m) == m: ok;  determinant(): ok
> checkMatrix(mc. <- as(cm, "Matrix"))  # dpo --> dsy -- (as above)  FIXME? ??
 norm(m [2 x 2]) : 1  I  F  M  ok
 Summary:  ok
2*m =?= m+m: identical
m >= m for all: ok
m < m for none: ok
symmpart(m) + skewpart(m) == m: ok;  determinant(): ok
> stopifnot(identical(mc., mcm),
+ 	  identical(.drop.factors(cm), (2*cm)/2),# remains dpo
+ 	  identical(cm + cp, cp + cs),# dge
+ 	  identical(mc., mcm),
+ 	  all(2*cm == mcm * 2))
> 
> checkMatrix(eq <- cm == cs)
 norm(m [2 x 2]) : 1  I  F  M  ok
 Summary:  ok
2*m =?= m+m: identical
m >= m for all: ok
m < m for none: ok
symmpart(m) + skewpart(m) == m: ok;  determinant(): ok
> stopifnot(all(eq@x),
+ 	  identical3(pack(eq), cs == cp, cm == cp),
+ 	  as.logical(!(cs < cp)),
+ 	  identical4(!(cs < cp), !(cp > cs), cp <= cs, cs >= cp))
> 
> ## Coercion to 'dpo' should give an error if result would be invalid
> M <- Matrix(diag(4) - 1)
> assertError(as(M, "dpoMatrix"))
> M. <- as(M, "generalMatrix")
> M.[1,2] <- 10 # -> not even symmetric anymore
> assertError(as(M., "dpoMatrix"))
> 
> 
> ## Cholesky
> checkMatrix(ch <- chol(cm))
 norm(m [2 x 2]) : 1  I  F  M  ok
 Summary:  ok
2*m =?= m+m: identical
m >= m for all: ok
m < m for none: ok
symmpart(m) + skewpart(m) == m: ok;  determinant(): ok
as(mm., "triangularMatrix"): valid: TRUE 
> checkMatrix(ch2 <- chol(as(cm, "dsyMatrix")))
 norm(m [2 x 2]) : 1  I  F  M  ok
 Summary:  ok
2*m =?= m+m: identical
m >= m for all: ok
m < m for none: ok
symmpart(m) + skewpart(m) == m: ok;  determinant(): ok
as(mm., "triangularMatrix"): valid: TRUE 
> checkMatrix(ch3 <- chol(as(cm, "generalMatrix")))
 norm(m [2 x 2]) : 1  I  F  M  ok
 Summary:  ok
2*m =?= m+m: identical
m >= m for all: ok
m < m for none: ok
symmpart(m) + skewpart(m) == m: ok;  determinant(): ok
as(mm., "triangularMatrix"): valid: TRUE 
> stopifnot(is.all.equal3(as(ch, "matrix"), as(ch2, "matrix"), as(ch3, "matrix")))
> ### Very basic	triangular matrix stuff
> 
> assertError( new("dtrMatrix", Dim = c(2,2), x= 1:4) )# double 'Dim'
> assertError( new("dtrMatrix", Dim = as.integer(c(2,2)), x= 1:4) )# int 'x'
> ## This caused a segfault (before revision r1172 in ../src/dtrMatrix.c):
> assertError( new("dtrMatrix", Dim = 2:2, x=as.double(1:4)) )# length(Dim) !=2
> assertError( new("dtrMatrix", Dim = as.integer(c(2,2)), x= as.double(1:5)))
> 
> tr22 <- new("dtrMatrix", Dim = as.integer(c(2,2)), x=as.double(1:4))
> tt22 <- t(tr22)
> (tPt <- tr22 + tt22)
2 x 2 Matrix of class "dgeMatrix"
     [,1] [,2]
[1,]    2    3
[2,]    3    8
> stopifnot(identical(10 * tPt, tPt * 10),
+ 	  as.vector(t.22 <- (tr22 / .5)* .5) == c(1,0,3,4),
+ 	  TRUE) ## not yet: class(t.22) == "dtrMatrix")
> 
> ## non-square triagonal Matrices --- are forbidden ---
> assertError(new("dtrMatrix", Dim = 2:3,
+ 		x=as.double(1:6), uplo="L", diag="U"))
> 
> n <- 3:3
> assertError(new("dtCMatrix", Dim = c(n,n), diag = "U"))
> validObject(T <- new("dtTMatrix", Dim = c(n,n), diag = "U"))
[1] TRUE
> validObject(M <- new("dtCMatrix", Dim = c(n,n), diag = "U",
+ 		     p = rep.int(0:0, n+1)))
[1] TRUE
> stopifnot(identical(as.mat(T), diag(n)))
> 
> suppressWarnings(RNGversion("3.5.0")); set.seed(3)
> (p9 <- p9. <- as(sample(9), "pMatrix"))
9 x 9 sparse Matrix of class "pMatrix"
                       
 [1,] . | . . . . . . .
 [2,] . . . . . . | . .
 [3,] . . | . . . . . .
 [4,] . . . . . . . . |
 [5,] . . . | . . . . .
 [6,] . . . . . . . | .
 [7,] | . . . . . . . .
 [8,] . . . . | . . . .
 [9,] . . . . . | . . .
> (np9 <- np9. <- as(p9, "TsparseMatrix"))
9 x 9 sparse Matrix of class "ngTMatrix"
                       
 [1,] . | . . . . . . .
 [2,] . . . . . . | . .
 [3,] . . | . . . . . .
 [4,] . . . . . . . . |
 [5,] . . . | . . . . .
 [6,] . . . . . . . | .
 [7,] | . . . . . . . .
 [8,] . . . . | . . . .
 [9,] . . . . . | . . .
> p9.[1,1] <- np9.[1,1] <- TRUE
> stopifnot(is.logical(p9[1,]),
+ 	  is(p9[2,, drop=FALSE], "indMatrix"),
+ 	  is(p9[9:1,], "indMatrix"),
+ 	  isTRUE(p9[-c(1:6, 8:9), 1]),
+ 	  identical(t(p9), solve(p9)),
+ 	  identical(p9[TRUE, ], p9),
+           identical(p9[, TRUE], p9),
+           identical(p9., np9.),
+ 	  identical(as(diag(9), "pMatrix"), as(1:9, "pMatrix"))
+ 	  )
> assert.EQ.mat(p9[TRUE,], as.matrix(np9))
> 
> ## validObject --> Cparse_validate(.)
> mm <- new("dgCMatrix", Dim = c(3L, 5L),
+           i = c(2L, 0L, 1L, 2L, 0L, 1L),
+           x = c( 2,  1,  1,  2,  1,  2),
+           p = c(0:2, 4L, 4L, 6L))
> 
> ## Previously unsorted columns were sorted - now are flagged as invalid
> m. <- mm
> ip <- c(1:2, 4:3, 6:5) # permute the 'i' and 'x' slot just "inside column":
> m.@i <- m.i <- mm@i[ip]
> m.@x <- m.x <- mm@x[ip]
> stopifnot(identical(1L, grep("not increasing within columns",
+                              validObject(m., test = TRUE))))
> .validateCsparse(m., TRUE) # don't use this at home, boys!
[1] TRUE
> m. # now is fixed
3 x 5 sparse Matrix of class "dgCMatrix"
              
[1,] . 1 . . 1
[2,] . . 1 . 2
[3,] 2 . 2 . .
> 
> ## Make sure that validObject() objects...
> ## 1) to wrong 'p'
> m. <- mm; m.@p[1] <- 1L
> stopifnot(identical(1L, grep("first element of 'p' slot",
+                              validObject(m., test = TRUE))))
> m.@p <- mm@p[c(1,3:2,4:6)]
> stopifnot(identical(1L, grep("not nondecreasing",
+                              validObject(m., test = TRUE))))
> ## 2) to non-strictly increasing i's:
> m. <- mm ; ix <- c(1:3,3,5:6)
> m.@i <- mm@i[ix]
> m.@x <- mm@x[ix]
> stopifnot(identical(1L, grep("not increasing within columns",
+                              validObject(m., test = TRUE))))
> ## ix <- c(1:3, 3:6) # now the the (i,x) slots are too large (and decreasing at end)
> ## m.@i <- mm@i[ix]
> ## m.@x <- mm@x[ix]
> ## stopifnot(identical(grep("^slot i is not.* increasing .*sort",
> ##                          (msg <- validObject(m., test=TRUE))),# seg.fault in the past
> ##                     1L))
> 
> ## over-allocation of the i- and x- slot should be allowed:
> ## (though it does not really help in M[.,.] <- *  yet)
> m. <- mm
> m.@i <- c(mm@i, NA, NA, NA)
> m.@x <- c(mm@x, 10:12)
> validObject(m.)
[1] TRUE
> m. # show() now works
3 x 5 sparse Matrix of class "dgCMatrix"
              
[1,] . 1 . . 1
[2,] . . 1 . 2
[3,] 2 . 2 . .
> stopifnot(all(m. == mm), # in spite of
+ 	  length(m.@i) > length(mm@i),
+           identical(t(t(m.)), mm),
+ 	  identical3(m. * m., m. * mm, mm * mm))
> m.[1,4] <- 99 ## FIXME: warning and cuts (!) the over-allocated slots
> 
> ## Low-level construction of invalid object:
> ##  Ensure that it does *NOT* segfault
> foo <- new("ngCMatrix",
+            i = as.integer(c(12204, 16799, 16799, 33517, 1128, 11930, 1128, 11930, 32183)),
+            p = rep(0:9, c(2,4,1,11,10,0,1,0,9,12)),
+            Dim = c(36952L, 49L))
> validObject(foo)# TRUE
[1] TRUE
> foo@i[5] <- foo@i[5] + 50000L
> msg <- validObject(foo, test=TRUE)# is -- correctly -- *not* valid anymore
> stopifnot(is.character(msg))
> ## Error in validObject(foo) :
> ##   invalid class "ngCMatrix" object: all row indices must be between 0 and nrow-1
> getLastMsg <- function(tryRes) {
+     ## Extract "final" message from erronous try result
+     sub("\n$", "", sub("^.*: ", "", as.character(tryRes)))
+ }
> t1 <- try(show(foo)) # error
36952 x 49 sparse Matrix of class "ngCMatrix"
Error in validObject(x) : 
  invalid class "ngCMatrix" object: 'i' slot has elements not in {0,...,Dim[1]-1}
> t2 <- try(head(foo)) # error, but previously a segfault
Error in ..subscript.2ary(x, l[[1L]], l[[2L]], drop = drop[1L]) : 
  invalid class "ngCMatrix" object: 'i' slot has elements not in {0,...,Dim[1]-1}
> stopifnot(identical(msg, getLastMsg(t1)),
+           identical(msg, getLastMsg(t2)))
> 
> ## test that all prototypes of nonvirtual classes are valid
> for (cl in c(names(getClassDef("Matrix")@subclasses),
+              names(getClassDef("MatrixFactorization")@subclasses)))
+     if(!isVirtualClass(def <- getClassDef(cl))) validObject(new(def))
> 
> cat('Time elapsed: ', proc.time(),'\n') # "stats"
Time elapsed:  3.14 0.26 3.39 NA NA 
> 
> if(!interactive()) warnings()
Warning message:
In Matrix(1:7 + 0, ncol = 3) :
  data length [7] is not a sub-multiple or multiple of the number of rows [3]
> 
> proc.time()
   user  system elapsed 
   3.14    0.26    3.39