# SCM Repository

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

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

revision 2136, Mon Mar 17 22:20:29 2008 UTC revision 2144, Tue Mar 18 23:08:12 2008 UTC
# Line 32  Line 32
32      }      }
33  }  }
34
35    ## Pkg 'spdep' had (relatively slow) versions of this as_dsCMatrix_I()
36    .symDiagonal <- function(n, x = rep.int(1,n), uplo = "U") {
37        stopifnot(n == (n. <- as.integer(n)), (n <- n.) >= 0)
38        if((lx <- length(x)) == 1) x <- rep.int(x, n)
39        else if(lx != n) stop("length(x) must be 1 or n")
40        cls <-
41            if(is.double(x)) "dsCMatrix"
42            else if(is.logical(x)) "lsCMatrix"
43            else { ## for now
44                storage.mode(x) <- "double"
45                "dsCMatrix"
46            }
47        new(cls, Dim = c(n,n), x = x, uplo = uplo,
48            i = if(n) 0:(n - 1L) else integer(0), p = 0:n)
49    }
50
51  ### This is modified from a post of Bert Gunter to R-help on  1 Sep 2005.  ### This is modified from a post of Bert Gunter to R-help on  1 Sep 2005.
52  ### Bert's code built on a post by Andy Liaw who most probably was influenced  ### Bert's code built on a post by Andy Liaw who most probably was influenced
53  ### by earlier posts, notably one by Scott Chasalow on S-news, 16 Jan 2002  ### by earlier posts, notably one by Scott Chasalow on S-news, 16 Jan 2002
# Line 88  Line 104
104  }  }
105
106  ## diagonal -> triangular,  upper / lower depending on "partner":  ## diagonal -> triangular,  upper / lower depending on "partner":
107  diag2tT.u <- function(d, x)  diag2tT.u <- function(d, x, kind = .M.kind(d))
108      .diag2tT(d, uplo = if(is(x,"triangularMatrix")) x@uplo else "U")      .diag2tT(d, uplo = if(is(x,"triangularMatrix")) x@uplo else "U", kind)
109
110
111  ## In order to evade method dispatch ambiguity warnings,  ## In order to evade method dispatch ambiguity warnings,
# Line 537  Line 553
553  ## For almost everything else, diag* shall be treated "as sparse" :  ## For almost everything else, diag* shall be treated "as sparse" :
554  ## These are cheap implementations via coercion  ## These are cheap implementations via coercion
555
556  ## for disambiguation --- define this for "sparseMatrix" , then for "ANY" :  ## For disambiguation --- define this for "sparseMatrix" , then for "ANY";
557  setMethod("Ops", signature(e1 = "diagonalMatrix", e2 = "sparseMatrix"),  ## and because we can save an .M.kind() call, we use this explicit
558            function(e1,e2) callGeneric(diag2tT.u(e1,e2), e2))  ## "hack" for all diagonalMatrix *subclasses* instead of just "diagonalMatrix" :
559  setMethod("Ops", signature(e1 = "sparseMatrix", e2 = "diagonalMatrix"),  ##
560            function(e1,e2) callGeneric(e1, diag2tT.u(e2,e1)))  ## ddi*:
561  ## in general:  setMethod("Ops", signature(e1 = "ddiMatrix", e2 = "sparseMatrix"),
562  setMethod("Ops", signature(e1 = "diagonalMatrix", e2 = "ANY"),            function(e1,e2) callGeneric(diag2tT.u(e1,e2, "d"), e2))
563            function(e1,e2) callGeneric(diag2tT.u(e1,e2), e2))  setMethod("Ops", signature(e1 = "sparseMatrix", e2 = "ddiMatrix"),
564  setMethod("Ops", signature(e1 = "ANY", e2 = "diagonalMatrix"),            function(e1,e2) callGeneric(e1, diag2tT.u(e2,e1, "d")))
565            function(e1,e2) callGeneric(e1, diag2tT.u(e2,e1)))  ## ldi*
566    setMethod("Ops", signature(e1 = "ldiMatrix", e2 = "sparseMatrix"),
567              function(e1,e2) callGeneric(diag2tT.u(e1,e2, "l"), e2))
568    setMethod("Ops", signature(e1 = "sparseMatrix", e2 = "ldiMatrix"),
569              function(e1,e2) callGeneric(e1, diag2tT.u(e2,e1, "l")))
570
571    ##  other = "numeric" : stay diagonal if possible
572    ## ddi*:
573    setMethod("Ops", signature(e1 = "ddiMatrix", e2 = "numeric"),
574              function(e1,e2) {
575                  n <- e1@Dim[1]; nsq <- n*n
576                  f0 <- callGeneric(0, e2)
577                  if(all(is0(f0))) { # remain diagonal
578                      if(e1@diag == "U" && (r <- callGeneric(1, e2)) != 1)
579                          e1@diag <- "N"
580                      else
581                          r <- callGeneric(e1@x, e2)
582                      e1@x <- if(length(e2) == nsq) r else rep(r, length.out = nsq)
583                      return(e1)
584                  }
585                  callGeneric(diag2tT.u(e1,e2, "d"), e2)
586              })
587
588    setMethod("Ops", signature(e1 = "numeric", e2 = "ddiMatrix"),
589              function(e1,e2) {
590                  n <- e2@Dim[1]; nsq <- n*n
591                  f0 <- callGeneric(e1, 0)
592                  if(all(is0(f0))) { # remain diagonal
593                      if(e2@diag == "U" && (r <- callGeneric(e1, 1)) != 1)
594                          e2@diag <- "N"
595                      else
596                          r <- callGeneric(e1, e2@x)
597                      e2@x <- if(length(e1) == nsq) r else rep(r, length.out = nsq)
598                      return(e2)
599                  }
600                  callGeneric(e1, diag2tT.u(e2,e1, "d"))
601              })
602    ## ldi*:
603    setMethod("Ops", signature(e1 = "ldiMatrix", e2 = "numeric"),
604              function(e1,e2) {
605                  n <- e1@Dim[1]; nsq <- n*n
606                  f0 <- callGeneric(FALSE, e2)
607                  if(all(is0(f0))) { # remain diagonal
608                      if(e1@diag == "U" && (r <- callGeneric(TRUE, e2)) != 1)
609                          e1@diag <- "N"
610                      else
611                          r <- callGeneric(e1@x, e2)
612                      e1@x <- if(length(e2) == nsq) r else rep(r, length.out = nsq)
613                      return(e1)
614                  }
615                  callGeneric(diag2tT.u(e1,e2, "l"), e2)
616              })
617
618    setMethod("Ops", signature(e1 = "numeric", e2 = "ldiMatrix"),
619              function(e1,e2) {
620                  n <- e2@Dim[1]; nsq <- n*n
621                  f0 <- callGeneric(e1, FALSE)
622                  if(all(is0(f0))) { # remain diagonal
623                      if(e2@diag == "U" && (r <- callGeneric(e1, TRUE)) != 1)
624                          e2@diag <- "N"
625                      else
626                          r <- callGeneric(e1, e2@x)
627                      e2@x <- if(length(e1) == nsq) r else rep(r, length.out = nsq)
628                      return(e2)
629                  }
630                  callGeneric(e1, diag2tT.u(e2,e1, "l"))
631              })
632
633
634
635    ## Not {"sparseMatrix", "numeric} :  {"denseMatrix", "matrix", ... }
636    ## ddi*:
637    setMethod("Ops", signature(e1 = "ddiMatrix", e2 = "ANY"),
638              function(e1,e2) callGeneric(diag2tT.u(e1,e2, "d"), e2))
639    setMethod("Ops", signature(e1 = "ANY", e2 = "ddiMatrix"),
640              function(e1,e2) callGeneric(e1, diag2tT.u(e2,e1, "d")))
641    ## ldi*:
642    setMethod("Ops", signature(e1 = "ldiMatrix", e2 = "ANY"),
643              function(e1,e2) callGeneric(diag2tT.u(e1,e2, "l"), e2))
644    setMethod("Ops", signature(e1 = "ANY", e2 = "ldiMatrix"),
645              function(e1,e2) callGeneric(e1, diag2tT.u(e2,e1, "l")))
646
647
648

Legend:
 Removed from v.2136 changed lines Added in v.2144