SCM

SCM Repository

[matrix] Diff of /pkg/R/Tsparse.R
ViewVC logotype

Diff of /pkg/R/Tsparse.R

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 1330, Fri Jul 21 08:28:18 2006 UTC revision 1331, Sat Jul 22 17:59:53 2006 UTC
# Line 1  Line 1 
1  #### "TsparseMatrix" : Virtual class of sparse matrices in triplet-format  #### "TsparseMatrix" : Virtual class of sparse matrices in triplet-format
2    
3  setAs("TsparseMatrix", "CsparseMatrix",  setAs("TsparseMatrix", "CsparseMatrix",
4        function(from) .Call(Tsparse_to_Csparse, from))        ## FIXME: this loses 'triangular' (and 'symmetric' ??)
5          function(from)
6          .Call(Tsparse_to_Csparse, from) ## ../src/Tsparse.c
7          ## |-> cholmod_T -> cholmod_C -> chm_sparse_to_SEXP(* , 1)
8          )
9    
10  ### "[" :  ### "[" :
11  ### -----  ### -----
# Line 146  Line 150 
150            if (drop && any(nd == 1)) drop(as(x,"matrix")) else x            if (drop && any(nd == 1)) drop(as(x,"matrix")) else x
151        })        })
152    
153    
154    ## workhorse for "[<-" :
155    replTmat <- function (x, i, j, value)
156    {
157        di <- dim(x)
158        dn <- dimnames(x)
159        i1 <- if(missing(i)) 0:(di[1] - 1:1) else .ind.prep2(i, 1, di, dn)
160        i2 <- if(missing(j)) 0:(di[2] - 1:1) else .ind.prep2(j, 2, di, dn)
161        dind <- c(length(i1), length(i2)) # dimension of replacement region
162        lenRepl <- prod(dind)
163        lenV <- length(value)
164        if(lenV == 0) {
165            if(lenRepl != 0)
166                stop("nothing to replace with")
167            else return(x)
168        }
169        ## else: lenV := length(value)       is > 0
170        if(lenRepl %% lenV != 0)
171            stop("number of items to replace is not a multiple of replacement length")
172    
173        ## Note: *T*matrix maybe non-unique: an entry can be split
174        ##    into a *sum* of several ones :
175        x <- uniq(x) # -> ./Auxiliaries.R
176    
177        sel <- ((m1 <- match(x@i, i1, nomatch=0)) > 0:0 &
178                (m2 <- match(x@j, i2, nomatch=0)) > 0:0)
179    
180        has.x <- any("x" == slotNames(x)) # i.e. *not* logical
181    
182        ## the simplest case: for all Tsparse, even for i or j missing
183        if(all(value == 0)) { ## just drop the non-zero entries
184            if(any(sel)) { ## non-zero there
185                x@i <- x@i[!sel]
186                x@j <- x@j[!sel]
187                if(has.x)
188                    x@x <- x@x[!sel]
189            }
190            return(x)
191        }
192    
193        ## else --  some( value != 0 ) --
194        if(lenV > lenRepl)
195            stop("too many replacement values")
196    
197        ## another simple, typical case:
198        if(lenRepl == 1) {
199            if(any(sel)) { ## non-zero there
200                if(has.x)
201                    x@x[sel] <- value
202            } else { ## new non-zero
203                x@i <- c(x@i, i1)
204                x@j <- c(x@j, i2)
205                if(has.x)
206                    x@x <- c(x@x, value)
207            }
208            return(x)
209        }
210    
211        v0 <- 0 == (value <- rep(value, length = lenRepl))
212        ## value[1:lenRepl]:  which are structural 0 now, which not?
213    
214        if(any(sel)) {
215            ## the 0-based indices of non-zero -- WRT to submatrix
216            non0 <- cbind(match(x@i[sel], i1),
217                          match(x@j[sel], i2)) - 1:1
218            iN0 <- 1:1 + encodeInd(non0, nr = dind[1])
219    
220            ## 1) replace those that are already non-zero (when value != 0)
221            vN0 <- !v0[iN0]
222            if(has.x)
223                x@x[sel[vN0]] <- value[iN0[vN0]]
224    
225            iI0 <- (1:lenRepl)[-iN0]        # == complementInd(non0, dind)
226        } else iI0 <- 1:lenRepl
227    
228        if(length(iI0)) {
229            ## 2) add those that were structural 0 (where value != 0)
230            vN0 <- !v0[iI0]
231            ij0 <- decodeInd(iI0[vN0] - 1:1, nr = dind[1])
232            x@i <- c(x@i, i1[ij0[,1] + 1:1])
233            x@j <- c(x@j, i2[ij0[,2] + 1:1])
234            if(has.x)
235                x@x <- c(x@x, value[iI0[vN0]])
236        }
237        x
238    }
239    
240    setReplaceMethod("[", signature(x = "TsparseMatrix", i = "index", j = "missing",
241                                    value = "replValue"),
242                     function (x, i, value) replTmat(x, i=i, value=value))
243    
244    setReplaceMethod("[", signature(x = "TsparseMatrix", i = "missing", j = "index",
245                                    value = "replValue"),
246                     function (x, j, value) replTmat(x, j=j, value=value))
247    
248    setReplaceMethod("[", signature(x = "TsparseMatrix", i = "index", j = "index",
249                                    value = "replValue"),
250                     replTmat)
251    
252    
253    
254    
255  setMethod("crossprod", signature(x = "TsparseMatrix", y = "missing"),  setMethod("crossprod", signature(x = "TsparseMatrix", y = "missing"),
256            function(x, y = NULL) {            function(x, y = NULL) {
257                a <- .Call(Csparse_crossprod, x, trans = FALSE, triplet = TRUE,                a <- .Call(Csparse_crossprod, x, trans = FALSE, triplet = TRUE)
                          PACKAGE = "Matrix")  
258                switch(substr(class(a)[1], 1, 1),                switch(substr(class(a)[1], 1, 1),
259                       "d" ={ new("dsCMatrix", i = a@i, p = a@p, x = a@x,                       "d" ={ new("dsCMatrix", i = a@i, p = a@p, x = a@x,
260                                  Dim = a@Dim, Dimnames = a@Dimnames, uplo = "U",                                  Dim = a@Dim, Dimnames = a@Dimnames, uplo = "U",
# Line 161  Line 266 
266    
267  setMethod("tcrossprod", signature(x = "TsparseMatrix", y = "missing"),  setMethod("tcrossprod", signature(x = "TsparseMatrix", y = "missing"),
268            function(x, y = NULL) {            function(x, y = NULL) {
269                a <- .Call(Csparse_crossprod, x, trans = TRUE, triplet = TRUE,                a <- .Call(Csparse_crossprod, x, trans = TRUE, triplet = TRUE)
                          PACKAGE = "Matrix")  
270                switch(substr(class(a)[1], 1, 1),                switch(substr(class(a)[1], 1, 1),
271                       "d" ={ new("dsCMatrix", i = a@i, p = a@p, x = a@x,                       "d" ={ new("dsCMatrix", i = a@i, p = a@p, x = a@x,
272                                  Dim = a@Dim, Dimnames = a@Dimnames, uplo = "L",                                  Dim = a@Dim, Dimnames = a@Dimnames, uplo = "L",

Legend:
Removed from v.1330  
changed lines
  Added in v.1331

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