SCM

SCM Repository

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

Annotation of /pkg/R/Rsparse.R

Parent Directory Parent Directory | Revision Log Revision Log


Revision 2207 - (view) (download)

1 : bates 677 #### Sparse Matrices in Compressed row-oriented format
2 : maechler 1599 #### --- "R"
3 : bates 677
4 : maechler 1734 ### ``mainly for completeness'' --- we *do* favour Csparse
5 :     ## - - - - - - - - - - - - hence only "minimal" methods here !
6 :     ## see also ./SparseM-conv.R
7 :    
8 : bates 677 ### contains = "dMatrix"
9 :    
10 : maechler 1760 ## compressed_to_TMatrix -- fails on 32bit--enable-R-shlib with segfault {Kurt}
11 :     ## ------------ --> ../src/dgCMatrix.c
12 : maechler 1766 .R.2.T <- function(from) .Call(compressed_to_TMatrix, from, FALSE)
13 : maechler 1760 ## slow R-level workaround
14 :     ## this is cheap; alternative: going there directly, using
15 :     ## i <- .Call(Matrix_expand_pointers, from@p),
16 : maechler 1766 if(FALSE)
17 : maechler 1760 .R.2.T <- function(from) as(.R.2.C(from), "TsparseMatrix")
18 : maechler 1758
19 : maechler 1911 ## R_to_CMatrix
20 : maechler 1758 ## ------------ --> ../src/dgCMatrix.c
21 : maechler 1766 .R.2.C <- function(from) .Call(R_to_CMatrix, from)
22 : maechler 1911
23 : maechler 2113 if(FALSE)## "slow" unneeded R-level version
24 : maechler 1758 .R.2.C <- function(from)
25 :     {
26 :     cl <- class(from)
27 :     valid <- c("dgRMatrix", "dsRMatrix", "dtRMatrix",
28 :     "lgRMatrix", "lsRMatrix", "ltRMatrix",
29 :     "ngRMatrix", "nsRMatrix", "ntRMatrix",
30 :     "zgRMatrix", "zsRMatrix", "ztRMatrix")
31 : maechler 1832 icl <- match(cl, valid) - 1L
32 : maechler 1758 if(is.na(icl)) stop("invalid class:", cl)
33 :     Ccl <- sub("^(..)R","\\1C", cl) # corresponding Csparse class name
34 :     r <- new(Ccl)
35 : maechler 1945 r@Dim <- from@Dim[2:1]
36 : maechler 1758 if(icl %/% 3 != 2) ## not "n..Matrix" --> has 'x' slot
37 :     r@x <- from@x
38 :     if(icl %% 3 != 0) { # symmetric or triangular
39 :     r@uplo <- from@uplo
40 :     if(icl %% 3 == 2) # triangular
41 :     r@diag <- from@diag
42 :     }
43 :     r@i <- from@j
44 :     r@p <- from@p
45 :     r <- t(r)
46 :     r@Dimnames <- from@Dimnames
47 :     r
48 :     }
49 : bates 677
50 : maechler 1911 ## However, a quick way to "treat a t(<R..>) as corresponding <C..> " :
51 :     .tR.2.C <- function(from)
52 :     {
53 :     cl <- class(from)
54 :     valid <- c("dgRMatrix", "dsRMatrix", "dtRMatrix",
55 :     "lgRMatrix", "lsRMatrix", "ltRMatrix",
56 :     "ngRMatrix", "nsRMatrix", "ntRMatrix",
57 :     "zgRMatrix", "zsRMatrix", "ztRMatrix")
58 :     icl <- match(cl, valid) - 1L
59 :     if(is.na(icl)) stop("invalid class:", cl)
60 :     Ccl <- sub("^(..)R","\\1C", cl) # corresponding Csparse class name
61 :     r <- new(Ccl)
62 :     r@i <- from@j
63 :     ##- -
64 :     r@p <- from@p
65 : maechler 1945 r@Dim <- from@Dim[2:1]
66 :     r@Dimnames <- from@Dimnames[2:1]
67 : maechler 1911
68 :     if(icl %/% 3 != 2) ## not "n..Matrix" --> has 'x' slot
69 :     r@x <- from@x
70 :     if(icl %% 3 != 0) { # symmetric or triangular
71 :     r@uplo <- from@uplo
72 :     if(icl %% 3 == 2) # triangular
73 :     r@diag <- from@diag
74 :     }
75 :     r
76 :     }
77 :    
78 :    
79 :    
80 : maechler 1760 ## coercion to other virtual classes --- the functionality we want to encourage
81 :    
82 : maechler 1751 setAs("RsparseMatrix", "TsparseMatrix", .R.2.T)
83 :     setAs("RsparseMatrix", "CsparseMatrix", .R.2.C)
84 : maechler 1329
85 : maechler 1760 setAs("RsparseMatrix", "denseMatrix",
86 :     function(from) as(.R.2.C(from), "denseMatrix"))
87 :    
88 : maechler 1751 setAs("RsparseMatrix", "dsparseMatrix",
89 : maechler 1760 function(from) as(.R.2.C(from), "dsparseMatrix"))
90 : maechler 1751 setAs("RsparseMatrix", "lsparseMatrix",
91 : maechler 1760 function(from) as(.R.2.C(from), "lsparseMatrix"))
92 : maechler 1751 setAs("RsparseMatrix", "nsparseMatrix",
93 : maechler 1760 function(from) as(.R.2.C(from), "nsparseMatrix"))
94 : maechler 1751
95 :     setAs("RsparseMatrix", "dMatrix",
96 : maechler 1760 function(from) as(.R.2.C(from), "dMatrix"))
97 : maechler 1751 setAs("RsparseMatrix", "lMatrix",
98 : maechler 1760 function(from) as(.R.2.C(from), "lMatrix"))
99 : maechler 1751 setAs("RsparseMatrix", "nMatrix",
100 : maechler 1760 function(from) as(.R.2.C(from), "nMatrix"))
101 : maechler 1751
102 : mmaechler 2175 setAs("RsparseMatrix", "generalMatrix",
103 :     function(from) as(.R.2.C(from), "generalMatrix"))
104 : maechler 1760
105 : mmaechler 2175
106 : maechler 1760 ## for printing etc:
107 :     setAs("RsparseMatrix", "dgeMatrix",
108 :     function(from) as(.R.2.C(from), "dgeMatrix"))
109 :     setAs("RsparseMatrix", "matrix",
110 :     function(from) as(.R.2.C(from), "matrix"))
111 :    
112 : maechler 2113 ## **VERY** cheap substitute: work via dgC and t(.)
113 : maechler 1751 .viaC.to.dgR <- function(from) {
114 : maechler 1332 m <- as(t(from), "dgCMatrix")
115 :     new("dgRMatrix", Dim = dim(from), Dimnames = .M.DN(from),
116 : maechler 1734 p = m@p, j = m@i, x = m@x)
117 : maechler 1332 }
118 :    
119 : maechler 2113 ## one of the few coercions "to <specific>" {tested in ../tests/Class+Meth.R}
120 :     setAs("matrix", "dgRMatrix", .viaC.to.dgR)
121 : maechler 1332
122 : maechler 2113 ## *very* cheap substitute: work via t(.) and Csparse
123 :     .viaC.to.R <- function(from) {
124 :     m <- as(t(from), "CsparseMatrix")# preserve symmetry/triangular
125 :     clx <- getClassDef(class(m))
126 :     has.x <- !extends(clx, "nsparseMatrix")## <==> has 'x' slot
127 :     ## instead of "d": .M.kind (m,cl)
128 :     ## instead of "g": ..M.shape(m,cl)
129 :     sh <- .M.shapeC(m,clx)
130 :     r <- new(paste(.M.kindC(clx), sh, "RMatrix", sep=""))
131 :     r@Dim <- dim(from)
132 :     r@Dimnames <- .M.DN(from)
133 :     r@p <- m@p
134 :     r@j <- m@i
135 :     if(has.x)
136 :     r@x <- m@x
137 :     if(sh != "g") {
138 : mmaechler 2175 r@uplo <- if(m@uplo != "U") "U" else "L"
139 : maechler 2113 if(sh == "t")
140 :     r@diag <- m@diag
141 :     }
142 :     r
143 :     }
144 :    
145 : mmaechler 2183 setAs("matrix", "RsparseMatrix", .viaC.to.R)
146 :     setAs("denseMatrix", "RsparseMatrix", .viaC.to.R)
147 :     setAs("sparseMatrix","RsparseMatrix", .viaC.to.R)
148 : maechler 2113
149 : maechler 1751 ## symmetric: can use same 'p' slot
150 : maechler 1734 setAs("dsCMatrix", "dsRMatrix",
151 :     function(from) new("dsRMatrix", Dim = dim(from), Dimnames = .M.DN(from),
152 :     p = from@p, j = from@i, x = from@x,
153 :     uplo = if (from@uplo == "U") "L" else "U"))
154 : maechler 1760 ## FIXME: if this makes sense, do it for "l" and "n" as well as "d"
155 : maechler 1332
156 : maechler 1760 ## setAs("dtCMatrix", "dtRMatrix", .viaC.to.dgR) # should work; can NOT use 'p'
157 : maechler 1332
158 : maechler 1734
159 : maechler 1332 ##setAs("dgRMatrix", "dgeMatrix",
160 :     ## function(from) .Call(csc_to_dgeMatrix, from))
161 :    
162 :     ##setAs("matrix", "dgRMatrix",
163 :     ## function(from) {
164 :     ## storage.mode(from) <- "double"
165 :     ## .Call(matrix_to_csc, from)
166 :     ## })
167 :    
168 :    
169 : maechler 1747 ##setMethod("diag", signature(x = "dgRMatrix"),
170 :     ## function(x = 1, nrow, ncol = n) .Call(csc_getDiag, x))
171 :    
172 : maechler 1189 ## try to define for "Matrix" -- once and for all -- but that fails -- why? __ FIXME __
173 :     ## setMethod("dim", signature(x = "dgRMatrix"),
174 :     ## function(x) x@Dim, valueClass = "integer")
175 : bates 677
176 :     ##setMethod("t", signature(x = "dgRMatrix"),
177 : maechler 1280 ## function(x) .Call(csc_transpose, x),
178 : bates 677 ## valueClass = "dgRMatrix")
179 :    
180 :     setMethod("image", "dgRMatrix",
181 :     function(x, ...) {
182 : maechler 1760 x <- as(x, "TsparseMatrix")
183 : bates 677 callGeneric()
184 :     })
185 : maechler 1349
186 : maechler 1760 setMethod("t", "RsparseMatrix", function(x) as(t(.R.2.T(x)), "RsparseMatrix"))
187 : maechler 1349
188 : maechler 1655
189 : maechler 1349 ## Want tril(), triu(), band() --- just as "indexing" ---
190 :     ## return a "close" class:
191 :     setMethod("tril", "RsparseMatrix",
192 : maechler 1760 function(x, k = 0, ...)
193 :     as(tril(.R.2.C(x), k = k, ...), "RsparseMatrix"))
194 : maechler 1349 setMethod("triu", "RsparseMatrix",
195 : maechler 1760 function(x, k = 0, ...)
196 :     as(triu(.R.2.C(x), k = k, ...), "RsparseMatrix"))
197 : maechler 1349 setMethod("band", "RsparseMatrix",
198 :     function(x, k1, k2, ...)
199 : maechler 1760 as(band(.R.2.C(x), k1 = k1, k2 = k2, ...), "RsparseMatrix"))
200 : maechler 2096
201 :     setReplaceMethod("[", signature(x = "RsparseMatrix", i = "index", j = "missing",
202 :     value = "replValue"),
203 :     function (x, i, j, ..., value)
204 :     replTmat(as(x,"TsparseMatrix"), i=i, value=value))
205 :    
206 :     setReplaceMethod("[", signature(x = "RsparseMatrix", i = "missing", j = "index",
207 :     value = "replValue"),
208 :     function (x, i, j, ..., value)
209 :     replTmat(as(x,"TsparseMatrix"), j=j, value=value))
210 :    
211 :     setReplaceMethod("[", signature(x = "RsparseMatrix", i = "index", j = "index",
212 :     value = "replValue"),
213 :     function (x, i, j, ..., value)
214 :     replTmat(as(x,"TsparseMatrix"), i=i, j=j, value=value))
215 :    
216 : mmaechler 2207 setReplaceMethod("[", signature(x = "RsparseMatrix", i = "index", j = "missing",
217 :     value = "sparseVector"),
218 :     function (x, i, j, ..., value)
219 :     replTmat(as(x,"TsparseMatrix"), i=i, value=value))
220 :    
221 :     setReplaceMethod("[", signature(x = "RsparseMatrix", i = "missing", j = "index",
222 :     value = "sparseVector"),
223 :     function (x, i, j, ..., value)
224 :     replTmat(as(x,"TsparseMatrix"), j=j, value=value))
225 :    
226 :     setReplaceMethod("[", signature(x = "RsparseMatrix", i = "index", j = "index",
227 :     value = "sparseVector"),
228 :     function (x, i, j, ..., value)
229 :     replTmat(as(x,"TsparseMatrix"), i=i, j=j, value=value))
230 :    
231 :    
232 : maechler 2096 setReplaceMethod("[", signature(x = "RsparseMatrix", i = "matrix", j = "missing",
233 :     value = "replValue"),
234 :     function (x, i, j, ..., value)
235 :     .TM.repl.i.2col(as(x,"TsparseMatrix"), i=i, value=value))
236 :    
237 :    
238 :    

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