fastMisc.Rd“Semi-API” functions used internally by Matrix,
often to bypass S4 dispatch and avoid the associated overhead.
These are exported to provide this capability to expert users.
Typical users should continue to rely on S4 generic functions
to dispatch suitable methods, by calling,
e.g., as(., <class>) for coercions.
.M2kind(from, kind = ".", sparse = NA)
.M2gen(from, kind = ".")
.M2sym(from, ...)
.M2tri(from, ...)
.M2diag(from)
.M2v(from)
.M2m(from)
.M2unpacked(from)
.M2packed(from)
.M2C(from)
.M2R(from)
.M2T(from)
.M2V(from)
.m2V(from, kind = ".")
.sparse2dense(from, packed = FALSE)
.diag2dense(from, kind = ".", shape = "t", packed = FALSE, uplo = "U")
.ind2dense(from, kind = "n")
.m2dense(from, class = ".ge", uplo = "U", diag = "N", trans = FALSE)
.dense2sparse(from, repr = "C")
.diag2sparse(from, kind = ".", shape = "t", repr = "C", uplo = "U")
.ind2sparse(from, kind = "n", repr = ".")
.m2sparse(from, class = ".gC", uplo = "U", diag = "N", trans = FALSE)
.tCRT(x, lazy = TRUE)
.diag.dsC(x, Chx = Cholesky(x, LDL = TRUE), res.kind = "diag")
.solve.dgC.lu (a, b, tol = .Machine$double.eps, check = TRUE)
.solve.dgC.qr (a, b, order = 3L, check = TRUE)
.solve.dgC.chol(a, b, check = TRUE)
.updateCHMfactor(object, parent, mult = 0)a Matrix, matrix, or vector.
a string (".", ",", "n", "l",
or "d") specifying the “kind” of the result.
"." indicates that the kind of from should be preserved.
"," is equivalent to "z" if from is complex
and to "d" otherwise.
"n" indicates that the result should inherit from
nMatrix or nsparseVector
(and so on).
a string (".", "g", "s", or
"t") specifying the “shape” of the result. "."
indicates that the shape of from should be preserved.
"g" indicates that the result should inherit from
generalMatrix (and so on).
a string (".", "C", "R", or
"T") specifying the sparse representation of the result.
"." is accepted only by .ind2sparse and indicates
the most efficient representation,
which is "C" ("R") for margin = 2 (1).
"C" indicates that the result should inherit from
CsparseMatrix (and so on).
a logical indicating if the result should
inherit from packedMatrix
rather than from unpackedMatrix.
It is ignored for from inheriting from
generalMatrix.
a logical indicating if the result should inherit
from sparseMatrix rather than from
denseMatrix. If NA, then the result
will be formally sparse if and only if from is.
a string ("U" or "L") indicating whether
the result (if symmetric or triangular) should store the upper or
lower triangle of from. The elements of from in the
opposite triangle are ignored.
a string ("N" or "U") indicating whether
the result (if triangular) should be formally nonunit or unit
triangular. In the unit triangular case, the diagonal elements
of from are ignored.
a logical indicating if the result should be a 1-row
matrix rather than a 1-column matrix where from is a vector
but not a matrix.
a string whose first three characters specify the class
of the result. It should match the pattern
"^[.nld](ge|sy|tr|sp|tp)" for .m2dense and
"^[.nld][gst][CRT]" for .m2sparse,
where "." in the first position is equivalent to "l"
for logical arguments and "d" for numeric arguments.
optional arguments passed to isSymmetric
or isTriangular.
a logical indicating if the transpose should be constructed with minimal allocation, but possibly without preserving representation.
optionally, the Cholesky(x, ...)
factorization of x. If supplied, then x is unused.
a string in c("trace", "sumLog", "prod", "min",
"max", "range", "diag", "diagBack").
see lu-methods.
see qr-methods.
a logical indicating if the first argument should be
tested for inheritance from dgCMatrix and
coerced if necessary. Set to FALSE for speed only if it
is known to already inherit from dgCMatrix.
Functions with names of the form .<A>2<B> implement coercions
from virtual class A to the “nearest” non-virtual subclass of
virtual class B, where the virtual classes are abbreviated as follows:
MVmmatrix
vvector
denseunpackedpackedsparseCRTgensymtridiagindAbbreviations should be seen as a guide, rather than as an
exact description of behaviour. Notably, .m2dense,
.m2sparse, and .m2V accept vectors that are
not matrices.
.tCRT(x)If lazy = TRUE, then .tCRT constructs the transpose
of x using the most efficient representation,
which for CRT is RCT. If lazy = FALSE,
then .tCRT preserves the representation of x,
behaving as the corresponding methods for generic function t.
.diag.dsC(x).diag.dsC computes (or uses if Chx is supplied)
the Cholesky factorization of x as \(L D L'\) in order
to calculate one of several possible statistics from the diagonal
entries of \(D\). See res.kind under ‘Arguments’.
.solve.dgC.*(a, b).solve.dgC.lu(a, b) needs a square matrix a.
.solve.dgC.qr(a, b) needs a “long” matrix a,
with nrow(a) >= ncol(a).
.solve.dgC.chol(a, b) needs a “wide” matrix a,
with nrow(a) <= ncol(a).
All three may be used to solve sparse linear systems directly.
Only .solve.dgC.qr and .solve.dgC.chol be used
to solve sparse least squares problems.
.updateCHMfactor(object, parent, mult).updateCHMfactor updates object with the result
of Cholesky factorizing
F(parent) + mult[1] * diag(nrow(parent)),
i.e., F(parent) plus mult[1] times the identity matrix,
where F = identity if parent is a dsCMatrix
and F = tcrossprod if parent is a dgCMatrix.
The nonzero pattern of F(parent) must match
that of S if object = Cholesky(S, ...).
D. <- diag(x = c(1, 1, 2, 3, 5, 8))
D.0 <- Diagonal(x = c(0, 0, 0, 3, 5, 8))
S. <- toeplitz(as.double(1:6))
C. <- new("dgCMatrix", Dim = c(3L, 4L),
p = c(0L, 1L, 1L, 1L, 3L), i = c(1L, 0L, 2L), x = c(-8, 2, 3))
stopifnot(exprs = {
identical(.M2tri (D.), as(D., "triangularMatrix"))
identical(.M2sym (D.), as(D., "symmetricMatrix"))
identical(.M2diag(D.), as(D., "diagonalMatrix"))
identical(.M2kind(C., "l"),
as(C., "lMatrix"))
identical(.M2kind(.sparse2dense(C.), "l"),
as(as(C., "denseMatrix"), "lMatrix"))
identical(.diag2sparse(D.0, ".", "t", "C"),
.dense2sparse(.diag2dense(D.0, ".", "t", TRUE), "C"))
identical(.M2gen(.diag2dense(D.0, ".", "s", FALSE)),
.sparse2dense(.M2gen(.diag2sparse(D.0, ".", "s", "T"))))
identical(S.,
.M2m(.m2sparse(S., ".sR")))
identical(S. * lower.tri(S.) + diag(1, 6L),
.M2m(.m2dense (S., ".tr", "L", "U")))
identical(.M2R(C.), .M2R(.M2T(C.)))
identical(.tCRT(C.), .M2R(t(C.)))
})
A <- tcrossprod(C.)/6 + Diagonal(3, 1/3); A[1,2] <- 3; A
#> 3 x 3 sparse Matrix of class "dgCMatrix"
#>
#> [1,] 1 3 1.000000
#> [2,] . 11 .
#> [3,] 1 . 1.833333
stopifnot(exprs = {
is.numeric( x. <- c(2.2, 0, -1.2) )
all.equal(x., .solve.dgC.lu(A, c(1,0,0), check=FALSE))
all.equal(x., .solve.dgC.qr(A, c(1,0,0), check=FALSE))
})
## Solving sparse least squares:
X <- rbind(A, Diagonal(3)) # design matrix X (for L.S.)
Xt <- t(X) # *transposed* X (for L.S.)
(y <- drop(crossprod(Xt, 1:3)) + c(-1,1)/1000) # small rand.err.
#> [1] 9.999 22.001 6.499 1.001 1.999 3.001
str(solveCh <- .solve.dgC.chol(Xt, y, check=FALSE)) # Xt *is* dgC..
#> List of 4
#> $ L :Formal class 'dCHMsimpl' [package "Matrix"] with 11 slots
#> .. ..@ x : num [1:6] 3 1 0.9444 128 0.0013 ...
#> .. ..@ p : int [1:4] 0 3 5 6
#> .. ..@ i : int [1:6] 0 1 2 1 2 2
#> .. ..@ nz : int [1:3] 3 2 1
#> .. ..@ nxt : int [1:5] 1 2 3 -1 0
#> .. ..@ prv : int [1:5] 4 0 1 2 -1
#> .. ..@ type : int [1:6] 2 0 0 1 0 0
#> .. ..@ colcount: int [1:3] 3 2 1
#> .. ..@ perm : int [1:3] 0 1 2
#> .. ..@ Dim : int [1:2] 3 3
#> .. ..@ Dimnames:List of 2
#> .. .. ..$ : NULL
#> .. .. ..$ : NULL
#> $ coef : num [1:3] 1 2 3
#> $ Xty : num [1:3] 17.5 274 24.9
#> $ resid: num [1:6] -0.000774 0.000308 -0.000306 0.00108 -0.001063 ...
stopifnot(exprs = {
all.equal(solveCh$coef, 1:3, tol = 1e-3)# rel.err ~ 1e-4
all.equal(solveCh$coef, drop(solve(tcrossprod(Xt), Xt %*% y)))
all.equal(solveCh$coef, .solve.dgC.qr(X, y, check=FALSE))
})