https://github.com/cran/Matrix
Tip revision: 1984a1439a1b9666a95fc64bc655504667e69e7b authored by Douglas Bates on 21 March 2006, 00:00:00 UTC
version 0.995-8
version 0.995-8
Tip revision: 1984a14
dtCMatrix.R
setMethod("t", signature(x = "dtCMatrix"),
function(x) {
tg <- t(as(x, "dgCMatrix"))
new("dtCMatrix", Dim = tg@Dim, Dimnames = tg@Dimnames,
p = tg@p, i = tg@i, x = tg@x, diag = x@diag,
uplo = ifelse(x@uplo == "U", "L", "U"))
}, valueClass = "dtCMatrix")
setAs("dtCMatrix", "ltCMatrix", # just drop 'x' slot:
function(from) new("ltCMatrix", i = from@i, p = from@p,
uplo = from@uplo, diag = from@diag,
## FIXME?: use from@factors smartly
Dim = from@Dim, Dimnames = from@Dimnames))
setAs("matrix", "dtCMatrix",
function(from) as(as(from, "dtTMatrix"), "dtCMatrix"))
setAs("dtCMatrix", "dgCMatrix",
function(from) {
if(from@diag == "U") { ## add diagonal of 1's
##FIXME: do this smartly - directly {in C or R}
as(as(from, "dgTMatrix"), "dgCMatrix")
}
else
new("dgCMatrix",
i = from@i, p = from@p, x = from@x,
Dim = from@Dim, Dimnames = from@Dimnames)
})
setAs("dtCMatrix", "dgTMatrix",
function(from)
.Call("tsc_to_dgTMatrix", from, PACKAGE = "Matrix"))
setAs("dtCMatrix", "dgeMatrix",
function(from) as(as(from, "dgTMatrix"), "dgeMatrix"))
## These are all needed because cholmod doesn't support triangular:
## (see end of ./Csparse.R )
setAs("dtCMatrix", "dtTMatrix",
function(from) {# and this is not elegant:
x <- as(from, "dgTMatrix")
## FIXME: if(from@diag == "U") should drop diagonal entries:
new("dtTMatrix", x = x@x, i = x@i, j = x@j,
Dim = x@Dim, Dimnames = x@Dimnames,
uplo = from@uplo, diag = "N")
})
setAs("dtCMatrix", "TsparseMatrix", function(from) as(from, "dtTMatrix"))
setAs("dtCMatrix", "dtrMatrix",
function(from) as(as(from, "dtTMatrix"), "dtrMatrix"))