# SCM Repository

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

# Diff of /pkg/Matrix/R/diagMatrix.R

revision 2778, Mon Apr 16 19:04:10 2012 UTC revision 2811, Mon Jul 16 15:50:31 2012 UTC
# Line 389  Line 389
389  ## FIXME: this now fails because the "denseMatrix" methods come first in dispatch  ## FIXME: this now fails because the "denseMatrix" methods come first in dispatch
390  ## Only(?) current bug:  x[i] <- value  is wrong when  i is *vector*  ## Only(?) current bug:  x[i] <- value  is wrong when  i is *vector*
391  replDiag <- function(x, i, j, ..., value) {  replDiag <- function(x, i, j, ..., value) {
392      x <- as(x, "TsparseMatrix")      x <- as(x, "CsparseMatrix")# was "Tsparse.." till 2012-07
393      if(missing(i))      if(missing(i))
394          x[, j] <- value          x[, j] <- value
395      else if(missing(j)) { ##  x[i , ] <- v  *OR*   x[i] <- v      else if(missing(j)) { ##  x[i , ] <- v  *OR*   x[i] <- v
# Line 429  Line 429
429                                   if(any(value != one | is.na(value))) {                                   if(any(value != one | is.na(value))) {
430                                       x@diag <- "N"                                       x@diag <- "N"
431                                       x@x <- rep.int(one, x@Dim[1])                                       x@x <- rep.int(one, x@Dim[1])
432                                   }                                   } else return(x)
433                               }                               }
434                               x@x[ii] <- value                               x@x[ii] <- value
435                               x                               x
# Line 726  Line 726
726  ###---------------- <Ops> (<Arith>, <Logic>, <Compare> ) ----------------------  ###---------------- <Ops> (<Arith>, <Logic>, <Compare> ) ----------------------
727
728  ## Use function for several signatures, in order to evade  ## Use function for several signatures, in order to evade
## ambiguous dispatch for "ddi", since there's also Arith(ddense., ddense.)
729  diagOdiag <- function(e1,e2) {  diagOdiag <- function(e1,e2) {
730      ## result should also be diagonal _ if possible _      ## result should also be diagonal _ if possible _
731      r <- callGeneric(.diag.x(e1), .diag.x(e2)) # error if not "compatible"      r <- callGeneric(.diag.x(e1), .diag.x(e2)) # error if not "compatible"
# Line 757  Line 756
756          stopifnot(length(r) == n)          stopifnot(length(r) == n)
757          xx <- as.vector(matrix(rbind(r, matrix(r00,n,n)), n,n))          xx <- as.vector(matrix(rbind(r, matrix(r00,n,n)), n,n))
758          newcl <-          newcl <-
759              paste(if(isNum) "d" else if(isLog) {              paste0(if(isNum) "d" else if(isLog) {
760                  if(!any(is.na(r)) && !any(is.na(r00))) "n" else "l"                  if(!any(is.na(r)) && !any(is.na(r00))) "n" else "l"
761              } else stop("not yet implemented .. please report")              } else stop("not yet implemented .. please report"), "syMatrix")
,
"syMatrix", sep='')
762
763          new(newcl, Dim = e1@Dim, Dimnames = e1@Dimnames, x = xx)          new(newcl, Dim = e1@Dim, Dimnames = e1@Dimnames, x = xx)
764      }      }
# Line 783  Line 780
780  ## diagonal  o  symmetric   |-->  symmetric  ## diagonal  o  symmetric   |-->  symmetric
781  ##    {also when other is sparse: do these "here" --  ##    {also when other is sparse: do these "here" --
782  ##     before conversion to sparse, since that loses "diagonality"}  ##     before conversion to sparse, since that loses "diagonality"}
783    diagOtri <- function(e1,e2) {
784        ## result must be triangular
785        r <- callGeneric(d1 <- .diag.x(e1), diag(e2)) # error if not "compatible"
786        ## Check what happens with non-diagonals, i.e. (0 o 0), (FALSE o 0), ...:
787        e1.0 <- if(.n1 <- is.numeric(d1   )) 0 else FALSE
788        r00 <- callGeneric(e1.0, if(.n2 <- is.numeric(e2[0])) 0 else FALSE)
789        if(is0(r00)) { ##  r00 == 0 or FALSE --- result *is* triangular
790            diag(e2) <- r
791            ## check what happens "in the triangle"
792            e2.2 <- if(.n2) 2 else TRUE
793            if(!callGeneric(e1.0, e2.2) == e2.2) { # values "in triangle" can change:
794                n <- dim(e2)[1L]
795                it <- indTri(n, upper = (e2@uplo == "U"))
796                e2[it] <- callGeneric(e1.0, e2[it])
797            }
798            e2
799        }
800        else { ## result not triangular ---> general
801            rr <- as(e2, "generalMatrix")
802            diag(rr) <- r
803            rr
804        }
805    }
806
807
808    setMethod("Ops", signature(e1 = "diagonalMatrix", e2 = "triangularMatrix"),
809              diagOtri)
810    ## For the reverse,  Ops == "Arith" | "Compare" | "Logic"
811    ##   'Arith'  :=  '"+"', '"-"', '"*"', '"^"', '"%%"', '"%/%"', '"/"'
812    setMethod("Arith", signature(e1 = "triangularMatrix", e2 = "diagonalMatrix"),
813              function(e1,e2)
814          { ## this must only trigger for *dense* e1
815              switch(.Generic,
816                     "+" = .Call(dtrMatrix_addDiag, as(e1,"dtrMatrix"),   .diag.x(e2)),
817                     "-" = .Call(dtrMatrix_addDiag, as(e1,"dtrMatrix"), - .diag.x(e2)),
818                     "*" = {
819                         n <- e2@Dim[1L]
820                         d2 <- if(e2@diag == "U") { # unit-diagonal
821                             d <- rep.int(as1(e2@x), n)
822                             e2@x <- d
823                             e2@diag <- "N"
824                             d
825                         } else e2@x
826                         e2@x <- diag(e1) * d2
827                         e2
828                     },
829                     "^" = { ## will be dense ( as  <ANY> ^ 0 == 1 ):
830                         e1 ^ as(e2, "denseMatrix")
831                     },
832                     ## otherwise:
833                     callGeneric(e1, diag2Tsmart(e2,e1)))
834    })
835
836    ## Compare --> 'swap' (e.g.   e1 < e2   <==>  e2 > e1 ):
837    setMethod("Compare", signature(e1 = "triangularMatrix", e2 = "diagonalMatrix"),
838              .Cmp.swap)
839    ## '&' and "|'  are commutative:
840    setMethod("Logic", signature(e1 = "triangularMatrix", e2 = "diagonalMatrix"),
841              function(e1,e2) callGeneric(e2, e1))
842
843  ## For almost everything else, diag* shall be treated "as sparse" :  ## For almost everything else, diag* shall be treated "as sparse" :
844  ## These are cheap implementations via coercion  ## These are cheap implementations via coercion

Legend:
 Removed from v.2778 changed lines Added in v.2811