SCM

SCM Repository

[matrix] Annotation of /pkg/R/Csparse.R
ViewVC logotype

Annotation of /pkg/R/Csparse.R

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1331 - (view) (download)

1 : maechler 1331 #### Methods for the virtual class 'CsparseMatrix' of sparse matrices stored in
2 :     #### "column compressed" format.
3 :     #### -- many more specific things are e.g. in ./dgCMatrix.R
4 :    
5 :     ## FIXME: the is(*,"generalMatrix") test at least makes these work,
6 :     ## but they are still ``wrong'', since triangularity is lost :
7 :    
8 :     setAs("CsparseMatrix", "TsparseMatrix",
9 :     function(from) {
10 :     if(!is(from, "generalMatrix")) { ## e.g. for triangular | symmetric
11 :     if (is(from, "dMatrix")) from <- as(from, "dgCMatrix")
12 :     else if(is(from, "lMatrix")) from <- as(from, "lgCMatrix")
13 :     else if(is(from, "zMatrix")) from <- as(from, "zgCMatrix")
14 :     else stop("undefined method for class ", class(from))
15 :     }
16 :     .Call(Csparse_to_Tsparse, from) # ../src/Csparse.c
17 :     ## |-> cholmod_C -> cholmod_T -> chm_triplet_to_SEXP(* , 1)
18 :     })
19 :    
20 :     setAs("CsparseMatrix", "denseMatrix",
21 :     function(from) {
22 :     if(!is(from, "generalMatrix")) { ## e.g. for triangular | symmetric
23 :     if (is(from, "dMatrix")) from <- as(from, "dgCMatrix")
24 :     else if(is(from, "lMatrix")) from <- as(from, "lgCMatrix")
25 :     else if(is(from, "zMatrix")) from <- as(from, "zgCMatrix")
26 :     else stop("undefined method for class ", class(from))
27 :     }
28 :     .Call(Csparse_to_dense, from)
29 :     })
30 :    
31 :     ### Some group methods:
32 :    
33 :     setMethod("Arith",
34 :     signature(e1 = "CsparseMatrix", e2 = "numeric"),
35 :     function(e1, e2) {
36 :     if(length(e2) == 1) { ## e.g., Mat ^ a
37 :     f0 <- callGeneric(0, e2)
38 :     if(!is.na(f0) && f0 == 0.) { # remain sparse, symm., tri.,...
39 :     e1@x <- callGeneric(e1@x, e2)
40 :     return(e1)
41 :     }
42 :     }
43 :     ## all other (potentially non-sparse) cases: give up symm, tri,..
44 :     callGeneric(as(e1, paste(.M.kind(e1), "gCMatrix", sep='')), e2)
45 :     })
46 :    
47 :     ## The same, e1 <-> e2 :
48 :     setMethod("Arith",
49 :     signature(e1 = "numeric", e2 = "CsparseMatrix"),
50 :     function(e1, e2) {
51 :     if(length(e1) == 1) {
52 :     f0 <- callGeneric(e1, 0)
53 :     if(!is.na(f0) && f0 == 0.) {
54 :     e2@x <- callGeneric(e1, e2@x)
55 :     return(e2)
56 :     }
57 :     }
58 :     callGeneric(e1, as(e2, paste(.M.kind(e2), "gCMatrix", sep='')))
59 :     })
60 :    
61 :    
62 :     setMethod("Math",
63 :     signature(x = "CsparseMatrix"),
64 :     function(x) {
65 :     f0 <- callGeneric(0.)
66 :     if(!is.na(f0) && f0 == 0.) {
67 :     ## sparseness, symm., triang.,... preserved
68 :     x@x <- callGeneric(x@x)
69 :     x
70 :     } else { ## no sparseness
71 :     callGeneric(as_dense(x))
72 :     }
73 :     })
74 :    
75 :    
76 :    
77 :     ### workhorse for "[<-" -- both for d* and l* C-sparse matrices :
78 :     replCmat <- function (x, i, j, value)
79 :     {
80 :     di <- dim(x)
81 :     dn <- dimnames(x)
82 :     i1 <- if(missing(i)) 0:(di[1] - 1:1) else .ind.prep2(i, 1, di, dn)
83 :     i2 <- if(missing(j)) 0:(di[2] - 1:1) else .ind.prep2(j, 2, di, dn)
84 :     dind <- c(length(i1), length(i2)) # dimension of replacement region
85 :     lenRepl <- prod(dind)
86 :     lenV <- length(value)
87 :     if(lenV == 0) {
88 :     if(lenRepl != 0)
89 :     stop("nothing to replace with")
90 :     else return(x)
91 :     }
92 :     ## else: lenV := length(value) is > 0
93 :     if(lenRepl %% lenV != 0)
94 :     stop("number of items to replace is not a multiple of replacement length")
95 :     if(lenV > lenRepl)
96 :     stop("too many replacement values")
97 :    
98 :     clx <- c(class(x))
99 :     xj <- .Call(Matrix_expand_pointers, x@p)
100 :     sel <- (!is.na(match(x@i, i1)) &
101 :     !is.na(match( xj, i2)))
102 :    
103 :     has.x <- any("x" == slotNames(x)) # i.e. *not* logical
104 :    
105 :     if(sum(sel) == lenRepl) { ## all entries to be replaced are non-zero:
106 :     value <- rep(value, length = lenRepl)
107 :     ## Ideally we only replace them where value != 0 and drop the value==0
108 :     ## ones; but that would have to (?) go through dgT*
109 :     ## v0 <- 0 == value
110 :     ## if (lenRepl == 1) and v0 is TRUE, the following is not doing anything
111 :     ##- --> ./dgTMatrix.R and its replTmat()
112 :     ## x@x[sel[!v0]] <- value[!v0]
113 :     x@x[sel] <- value
114 :     return(x)
115 :     }
116 :     ## else go via Tsparse..
117 :     x <- as(x, "TsparseMatrix")
118 :     x[i,j] <- value
119 :     as_CspClass(x, clx)
120 :     }
121 :    
122 :     setReplaceMethod("[", signature(x = "CsparseMatrix", i = "index", j = "missing",
123 :     value = "replValue"),
124 :     function (x, i, value) replCmat(x, i=i, value=value))
125 :    
126 :     setReplaceMethod("[", signature(x = "CsparseMatrix", i = "missing", j = "index",
127 :     value = "replValue"),
128 :     function (x, j, value) replCmat(x, j=j, value=value))
129 :    
130 :     setReplaceMethod("[", signature(x = "CsparseMatrix", i = "index", j = "index",
131 :     value = "replValue"),
132 :     replCmat)
133 :    
134 :    
135 : bates 923 setMethod("crossprod", signature(x = "CsparseMatrix", y = "missing"),
136 : bates 1267 function(x, y = NULL) {
137 : maechler 1331 a <- .Call(Csparse_crossprod, x, trans = FALSE, triplet = FALSE)
138 : maechler 1270 switch(substr(class(a)[1], 1, 1),
139 :     "d" ={ new("dsCMatrix", i = a@i, p = a@p, x = a@x,
140 :     Dim = a@Dim, Dimnames = a@Dimnames, uplo = "U",
141 :     factors = list()) },
142 :     "l" ={ new("lsCMatrix", i = a@i, p = a@p,
143 :     Dim = a@Dim, Dimnames = a@Dimnames, uplo = "U",
144 :     factors = list()) })
145 :     })
146 : bates 923
147 : bates 1267
148 : bates 923 setMethod("t", signature(x = "CsparseMatrix"),
149 : maechler 958 function(x)
150 : maechler 1280 .Call(Csparse_transpose, x))
151 : maechler 958
152 : maechler 1087 setMethod("tcrossprod", signature(x = "CsparseMatrix", y = "missing"),
153 : bates 1267 function(x, y = NULL) {
154 : maechler 1331 a <- .Call(Csparse_crossprod, x, trans = TRUE, triplet = FALSE)
155 : maechler 1270 switch(substr(class(a)[1], 1, 1),
156 :     "d" ={ new("dsCMatrix", i = a@i, p = a@p, x = a@x,
157 :     Dim = a@Dim, Dimnames = a@Dimnames, uplo = "L",
158 :     factors = list()) },
159 :     "l" ={ new("lsCMatrix", i = a@i, p = a@p,
160 :     Dim = a@Dim, Dimnames = a@Dimnames, uplo = "L",
161 :     factors = list()) })
162 :     })
163 : bates 1267
164 : maechler 1087 ## FIXME (TODO):
165 :     ## setMethod("tcrossprod", signature(x = "CsparseMatrix", y = "CsparseMatrix"),
166 :     ## function(x, y)
167 : maechler 1331 ## .Call(Csparse_crossprod_2, x, y, trans = TRUE, triplet = FALSE)
168 : bates 1059
169 : maechler 1087
170 : bates 1059 setMethod("%*%", signature(x = "CsparseMatrix", y = "CsparseMatrix"),
171 : maechler 1280 function(x, y) .Call(Csparse_Csparse_prod, x, y))
172 : bates 1059
173 :     setMethod("%*%", signature(x = "CsparseMatrix", y = "denseMatrix"),
174 : maechler 1280 function(x, y) .Call(Csparse_dense_prod, x, y))
175 : bates 1059
176 : maechler 1331 ## NB: have extra tril(), triu() methods for "dsC" and "lsC"
177 : bates 1268 setMethod("tril", "CsparseMatrix",
178 : maechler 1331 function(x, k = 0, ...) {
179 :     k <- as.integer(k[1])
180 :     dd <- dim(x)
181 :     stopifnot(-dd[1] <= k, k <= dd[1]) # had k <= 0
182 :     r <- .Call(Csparse_band, x, -dd[1], k)
183 :     ## return "lower triangular" if k <= 0
184 :     if(k <= 0) as(r, paste(.M.kind(x), "tCMatrix", sep='')) else r
185 :     })
186 : bates 1265
187 : bates 1268 setMethod("triu", "CsparseMatrix",
188 : maechler 1331 function(x, k = 0, ...) {
189 :     k <- as.integer(k[1])
190 :     dd <- dim(x)
191 :     stopifnot(-dd[1] <= k, k <= dd[1]) # had k >= 0
192 :     r <- .Call(Csparse_band, x, k, dd[2])
193 :     ## return "upper triangular" if k >= 0
194 :     if(k >= 0) as(r, paste(.M.kind(x), "tCMatrix", sep='')) else r
195 :     })
196 : bates 1268
197 :     setMethod("band", "CsparseMatrix",
198 : maechler 1331 function(x, k1, k2, ...) {
199 :     k1 <- as.integer(k1[1])
200 :     k2 <- as.integer(k2[1])
201 :     dd <- dim(x)
202 :     stopifnot(-dd[1] <= k1, k1 <= k2, k2 <= dd[1])
203 :     r <- .Call(Csparse_band, x, k1, k2)
204 :     if(k1 * k2 >= 0) ## triangular
205 :     as(r, paste(.M.kind(x), "tCMatrix", sep=''))
206 :     else if (k1 < 0 && k1 == -k2 && isSymmetric(x)) ## symmetric
207 :     as(r, paste(.M.kind(x), "sCMatrix", sep=''))
208 :     else
209 :     r
210 :     })
211 : maechler 1290
212 :     setMethod("colSums", signature(x = "CsparseMatrix"), .as.dgC.Fun,
213 :     valueClass = "numeric")
214 :     setMethod("colMeans", signature(x = "CsparseMatrix"), .as.dgC.Fun,
215 :     valueClass = "numeric")
216 :     setMethod("rowSums", signature(x = "CsparseMatrix"), .as.dgC.Fun,
217 :     valueClass = "numeric")
218 :     setMethod("rowMeans", signature(x = "CsparseMatrix"), .as.dgC.Fun,
219 :     valueClass = "numeric")

root@r-forge.r-project.org
ViewVC Help
Powered by ViewVC 1.0.0  
Thanks to:
Vienna University of Economics and Business Powered By FusionForge