# SCM Repository

[matrix] Diff of /pkg/R/Matrix.R
 [matrix] / pkg / R / Matrix.R

# Diff of /pkg/R/Matrix.R

revision 2110, Sat Jan 26 20:59:26 2008 UTC revision 2112, Mon Feb 18 08:24:46 2008 UTC
# Line 49  Line 49
49
50  ## "base" has an isSymmetric() S3-generic since R 2.3.0  ## "base" has an isSymmetric() S3-generic since R 2.3.0
51  setMethod("isSymmetric", signature(object = "symmetricMatrix"),  setMethod("isSymmetric", signature(object = "symmetricMatrix"),
52            function(object,tol) TRUE)            function(object, ...) TRUE)
53  setMethod("isSymmetric", signature(object = "triangularMatrix"),  setMethod("isSymmetric", signature(object = "triangularMatrix"),
54            ## TRUE iff diagonal:            ## TRUE iff diagonal:
55            function(object,tol) isDiagonal(object))            function(object, ...) isDiagonal(object))

setMethod("isTriangular", signature(object = "triangularMatrix"),
function(object, ...) TRUE)
56
57  setMethod("isTriangular", signature(object = "matrix"), isTriMat)  setMethod("isTriangular", signature(object = "matrix"), isTriMat)
58
59  setMethod("isDiagonal", signature(object = "matrix"), .is.diagonal)  setMethod("isDiagonal", signature(object = "matrix"), .is.diagonal)
60
61    ## The "catch all" methods -- far from optimal:
62    setMethod("symmpart", signature(x = "Matrix"),
63              function(x) as((x + t(x))/2, "symmetricMatrix"))
64    setMethod("skewpart", signature(x = "Matrix"),
65              function(x) (x - t(x))/2)
66
67    ## FIXME: do this (similarly as for "ddense.." in C
68    setMethod("symmpart", signature(x = "matrix"), function(x) (x + t(x))/2)
69    setMethod("skewpart", signature(x = "matrix"), function(x) (x - t(x))/2)
70
71
72
73
74  setMethod("dim", signature(x = "Matrix"),  setMethod("dim", signature(x = "Matrix"),
# Line 376  Line 384
384
385  ## missing 'drop' --> 'drop = TRUE'  ## missing 'drop' --> 'drop = TRUE'
386  ##                     -----------  ##                     -----------
387  ## select rows  ## select rows __ or __ vector indexing:
388  setMethod("[", signature(x = "Matrix", i = "index", j = "missing",  setMethod("[", signature(x = "Matrix", i = "index", j = "missing",
389                           drop = "missing"),                           drop = "missing"),
390            function(x,i,j, ..., drop) {            function(x,i,j, ..., drop) {
391                if(nargs() == 2) { ## e.g. M[0] , M[TRUE],  M[1:2]                if(nargs() == 2) { ## e.g. M[0] , M[TRUE],  M[1:2]
392                    if(any(i) || prod(dim(x)) == 0)                    if(any(as.logical(i)) || prod(dim(x)) == 0)
393                          ## FIXME: for *large sparse*, use sparseVector !
394                        as.vector(x)[i]                        as.vector(x)[i]
395                    else ## save memory                    else ## save memory (for large sparse M):
396                        as.vector(x[1,1])[FALSE]                        as.vector(x[1,1])[FALSE]
397                } else {                } else {
398                    callGeneric(x, i=i, , drop=TRUE)                    callGeneric(x, i=i, , drop=TRUE)
# Line 412  Line 421
421      nA <- nargs()      nA <- nargs()
422      if(nA == 2) { ##  M [ M >= 7 ]      if(nA == 2) { ##  M [ M >= 7 ]
423          ## FIXME: when both 'x' and 'i' are sparse, this can be very inefficient          ## FIXME: when both 'x' and 'i' are sparse, this can be very inefficient
424          as(x, geClass(x))@x[as.vector(i)]          if(is(x, "sparseMatrix"))
425                message("<sparse>[ <logic> ] : .M.sub.i.logical() maybe inefficient")
426            toC <- geClass(x)
427            if(canCoerce(x, toC)) as(x, toC)@x[as.vector(i)]
428            else as(as(as(x, "generalMatrix"), "denseMatrix"), toC)@x[as.vector(i)]
429          ## -> error when lengths don't match          ## -> error when lengths don't match
430      } else if(nA == 3) { ##  M [ M[,1, drop=FALSE] >= 7, ]      } else if(nA == 3) { ##  M [ M[,1, drop=FALSE] >= 7, ]
431          stop("not-yet-implemented 'Matrix' subsetting") ## FIXME          stop("not-yet-implemented 'Matrix' subsetting") ## FIXME
# Line 428  Line 441
441            .M.sub.i.logical)            .M.sub.i.logical)
442
443
444    subset.ij <- function(x, ij) {
445        m <- nrow(ij)
446        if(m > 3) {
447            cld <- getClassDef(class(x))
448            sym.x <- extends(cld, "symmetricMatrix")
449            if(sym.x) {
450                W <- if(x@uplo == "U") # stored only [i,j] with i <= j
451                    ij[,1] > ij[,2] else ij[,1] < ij[,2]
452                if(any(W))
453                    ij[W,] <- ij[W, 2:1]
454            }
455            if(extends(cld, "sparseMatrix")) {
456                ## do something smarter:
457                nr <- nrow(x)
458                if(!extends(cld, "CsparseMatrix")) {
459                    x <- as(x, "CsparseMatrix") # simpler; our standard
460                    cld <- getClassDef(class(x))
461                }
462                tri.x <- extends(cld, "triangularMatrix")
463                if(tri.x) {
464                    ## need these for the 'x' slot in any case
465                    if (x@diag == "U") x <- .Call(Csparse_diagU2N, x)
466                    ## slightly more efficient than non0.i() or non0ind():
467                    ij.x <- .Call(compressed_non_0_ij, x, isC=TRUE)
468                } else { ## symmetric / general : for symmetric, only "existing"b
469                    ij.x <- non0.i(x, cld)
470                }
471
472                mi <- match(encodeInd(ij.x,   nr),
473                            encodeInd(ij -1L, nr), nomatch=0)
474                mmi <- mi != 0
475                ## Result:
476                ans <- vector(mode = .type.kind[.M.kindC(cld)], length = m)
477                ## those that are *not* zero:
478                ans[mi[mmi]] <-
479                    if(extends(cld, "nsparseMatrix")) TRUE else x@x[mmi]
480                ans
481
482            } else { ## non-sparse : dense
483                message("m[ <ij-matrix> ]: inefficiently indexing single elements")
484                i1 <- ij[,1]
485                i2 <- ij[,2]
486                ## very inefficient for large m -- FIXME! --
487                unlist(lapply(seq_len(m), function(j) x[i1[j], i2[j]]))
488            }
489        } else { # 1 <= m <= 3
490            i1 <- ij[,1]
491            i2 <- ij[,2]
492            unlist(lapply(seq_len(m), function(j) x[i1[j], i2[j]]))
493        }
494    }
495
496  ## A[ ij ]  where ij is (i,j) 2-column matrix -- but also when that is logical mat!  ## A[ ij ]  where ij is (i,j) 2-column matrix -- but also when that is logical mat!
497  .M.sub.i.2col <- function (x, i, j, ..., drop)  .M.sub.i.2col <- function (x, i, j, ..., drop)
498  {  {
# Line 443  Line 508
508          m <- nrow(i)          m <- nrow(i)
509          if(m == 0) return(vector(mode = .type.kind[.M.kind(x)]))          if(m == 0) return(vector(mode = .type.kind[.M.kind(x)]))
510          ## else          ## else
511          i1 <- i[,1]          subset.ij(x, i)
i2 <- i[,2]
## potentially inefficient -- FIXME --
unlist(lapply(seq_len(m), function(j) x[i1[j], i2[j]]))
512
513      } else stop("nargs() = ", nA,      } else stop("nargs() = ", nA,
514                  ".  Extraneous illegal arguments inside '[ .. ]' (i.2col)?")                  ".  Extraneous illegal arguments inside '[ .. ]' (i.2col)?")
# Line 502  Line 564
564          value <- rep(value, length = m)          value <- rep(value, length = m)
565          i1 <- i[,1]          i1 <- i[,1]
566          i2 <- i[,2]          i2 <- i[,2]
567            if(m > 2)
568                message("m[ <ij-matrix> ] <- v: inefficiently treating single elements")
569          ## inefficient -- FIXME -- (also loses "symmetry" unnecessarily)          ## inefficient -- FIXME -- (also loses "symmetry" unnecessarily)
570          for(k in seq_len(m))          for(k in seq_len(m))
571              x[i1[k], i2[k]] <- value[k]              x[i1[k], i2[k]] <- value[k]

Legend:
 Removed from v.2110 changed lines Added in v.2112