SCM

SCM Repository

[matrix] Annotation of /pkg/Matrix/R/HBMM.R
ViewVC logotype

Annotation of /pkg/Matrix/R/HBMM.R

Parent Directory Parent Directory | Revision Log Revision Log


Revision 3046 - (view) (download)

1 : bates 825 ## Utilities for the Harwell-Boeing and MatrixMarket formats
2 :    
3 : bates 857 readone <- function(ln, iwd, nper, conv)
4 : bates 846 {
5 : bates 904 ln <- gsub("D", "E", ln)
6 : bates 857 inds <- seq(0, by = iwd, length = nper + 1)
7 :     (conv)(substring(ln, 1 + inds[-length(inds)], inds[-1]))
8 : bates 846 }
9 :    
10 : bates 857 readmany <- function(conn, nlines, nvals, fmt, conv)
11 : bates 846 {
12 : bates 857 if (!grep("[[:digit:]]+[DEFGI][[:digit:]]+", fmt))
13 : maechler 1604 stop("Not a valid format")
14 : bates 857 Iind <- regexpr('[DEFGI]', fmt)
15 :     nper <- as.integer(substr(fmt, regexpr('[[:digit:]]+[DEFGI]', fmt), Iind - 1))
16 : maechler 1604 iwd <- as.integer(substr(fmt, Iind + 1, regexpr('[\\.\\)]', fmt) - 1))
17 : bates 857 rem <- nvals %% nper
18 :     full <- nvals %/% nper
19 :     ans <- vector("list", nvals %/% nper)
20 : maechler 1654 for (i in seq_len(full))
21 : maechler 1604 ans[[i]] <- readone(readLines(conn, 1, ok = FALSE),
22 :     iwd, nper, conv)
23 : bates 857 if (!rem) return(unlist(ans))
24 :     c(unlist(ans),
25 :     readone(readLines(conn, 1, ok = FALSE), iwd, rem, conv))
26 : bates 846 }
27 :    
28 : bates 852 readHB <- function(file)
29 : bates 846 {
30 : maechler 1604 if (is.character(file))
31 :     file <- if (file == "") stdin() else file(file)
32 :     if (!inherits(file, "connection"))
33 : bates 852 stop("'file' must be a character string or connection")
34 :     if (!isOpen(file)) {
35 :     open(file)
36 :     on.exit(close(file))
37 :     }
38 :     hdr <- readLines(file, 4, ok = FALSE)
39 : mmaechler 3046 ## Title <- sub('[[:space:]]+$', '', substr(hdr[1], 1, 72))
40 :     ## Key <- sub('[[:space:]]+$', '', substr(hdr[1], 73, 80))
41 :     ## totln <- as.integer(substr(hdr[2], 1, 14))
42 : bates 846 ptrln <- as.integer(substr(hdr[2], 15, 28))
43 :     indln <- as.integer(substr(hdr[2], 29, 42))
44 :     valln <- as.integer(substr(hdr[2], 43, 56))
45 :     rhsln <- as.integer(substr(hdr[2], 57, 70))
46 :     if (!(t1 <- substr(hdr[3], 1, 1)) %in% c('C', 'R', 'P'))
47 : mmaechler 2863 stop(gettextf("Invalid storage type: %s", t1), domain=NA)
48 : bates 846 if (t1 != 'R') stop("Only numeric sparse matrices allowed")
49 :     ## _FIXME: Patterns should also be allowed
50 :     if (!(t2 <- substr(hdr[3], 2, 2)) %in% c('H', 'R', 'S', 'U', 'Z'))
51 : mmaechler 2863 stop(gettextf("Invalid storage format: %s", t2), domain=NA)
52 : bates 846 if (!(t3 <- substr(hdr[3], 3, 3)) %in% c('A', 'E'))
53 : mmaechler 2863 stop(gettextf("Invalid assembled indicator: %s", t3), domain=NA)
54 : bates 846 nr <- as.integer(substr(hdr[3], 15, 28))
55 :     nc <- as.integer(substr(hdr[3], 29, 42))
56 :     nz <- as.integer(substr(hdr[3], 43, 56))
57 : mmaechler 3046 ## nel <- as.integer(substr(hdr[3], 57, 70))
58 : bates 852 ptrfmt <- toupper(sub('[[:space:]]+$', '', substr(hdr[4], 1, 16)))
59 :     indfmt <- toupper(sub('[[:space:]]+$', '', substr(hdr[4], 17, 32)))
60 :     valfmt <- toupper(sub('[[:space:]]+$', '', substr(hdr[4], 33, 52)))
61 : mmaechler 3046 ## rhsfmt <- toupper(sub('[[:space:]]+$', '', substr(hdr[4], 53, 72)))
62 :     if (!is.na(rhsln) && rhsln > 0) readLines(file, 1, ok = FALSE) # h5
63 : bates 857 ptr <- readmany(file, ptrln, nc + 1, ptrfmt, as.integer)
64 :     ind <- readmany(file, indln, nz, indfmt, as.integer)
65 :     vals <- readmany(file, valln, nz, valfmt, as.numeric)
66 : bates 846 if (t2 == 'S')
67 : maechler 1832 new("dsCMatrix", uplo = "L", p = ptr - 1L,
68 :     i = ind - 1L, x = vals, Dim = c(nr, nc))
69 : maechler 1604 else
70 : maechler 1832 new("dgCMatrix", p = ptr - 1L,
71 :     i = ind - 1L, x = vals, Dim = c(nr, nc))
72 : bates 825
73 : bates 846 }
74 :    
75 : bates 852 readMM <- function(file)
76 :     {
77 : maechler 1604 if (is.character(file))
78 : maechler 2054 file <- if(file == "") stdin() else file(file)
79 : maechler 1604 if (!inherits(file, "connection"))
80 : maechler 2054 stop("'file' must be a character string or connection")
81 : bates 852 if (!isOpen(file)) {
82 : maechler 2054 open(file)
83 :     on.exit(close(file))
84 : bates 852 }
85 : maechler 2054 scan1 <- function(what, ...)
86 :     scan(file, nmax = 1, what = what, quiet = TRUE, ...)
87 :    
88 : mmaechler 3046 if (scan1(character()) != "%%MatrixMarket")# hdr
89 : maechler 2054 stop("file is not a MatrixMarket file")
90 :     if (!(typ <- tolower(scan1(character()))) %in% "matrix")
91 : mmaechler 2822 stop(gettextf("type '%s' not recognized", typ), domain = NA)
92 : maechler 2054 if (!(repr <- tolower(scan1(character()))) %in% c("coordinate", "array"))
93 : mmaechler 2822 stop(gettextf("representation '%s' not recognized", repr), domain = NA)
94 : maechler 2054 elt <- tolower(scan1(character()))
95 : bates 852 if (!elt %in% c("real", "complex", "integer", "pattern"))
96 : mmaechler 2822 stop(gettextf("element type '%s' not recognized", elt), domain = NA)
97 :    
98 : maechler 2054 sym <- tolower(scan1(character()))
99 : bates 852 if (!sym %in% c("general", "symmetric", "skew-symmetric", "hermitian"))
100 : mmaechler 2822 stop(gettextf("symmetry form '%s' not recognized", sym), domain = NA)
101 : maechler 2054 nr <- scan1(integer(), comment.char = "%")
102 :     nc <- scan1(integer())
103 :     nz <- scan1(integer())
104 : mmaechler 2369 checkIJ <- function(els) {
105 :     if(els$i < 1 || els$i > nr)
106 :     stop("readMM(): row values 'i' are not in 1:nr", call.=FALSE)
107 :     if(els$j < 1 || els$j > nc)
108 :     stop("readMM(): column values 'j' are not in 1:nc", call.=FALSE)
109 :     }
110 : maechler 2054 if (repr == "coordinate") {
111 :     switch(elt,
112 :     "real" = ,
113 :     "integer" = {
114 :     ## TODO: the "integer" element type should be returned as
115 :     ## an object of an "iMatrix" subclass--once there are
116 :     els <- scan(file, nmax = nz, quiet = TRUE,
117 :     what= list(i= integer(), j= integer(), x= numeric()))
118 : mmaechler 2369 checkIJ(els)
119 : maechler 2054 switch(sym,
120 :     "general" = {
121 :     new("dgTMatrix", Dim = c(nr, nc), i = els$i - 1L,
122 :     j = els$j - 1L, x = els$x)
123 :     },
124 :     "symmetric" = {
125 :     new("dsTMatrix", uplo = "L", Dim = c(nr, nc),
126 :     i = els$i - 1L, j = els$j - 1L, x = els$x)
127 :     },
128 :     "skew-symmetric" = {
129 :     stop("symmetry form 'skew-symmetric' not yet implemented for reading")
130 :     ## FIXME: use dgT... but must expand the (i,j,x) slots!
131 :     new("dgTMatrix", uplo = "L", Dim = c(nr, nc),
132 :     i = els$i - 1L, j = els$j - 1L, x = els$x)
133 :    
134 :     },
135 :     "hermitian" = {
136 :     stop("symmetry form 'hermitian' not yet implemented for reading")
137 :     },
138 :     ## otherwise (not possible; just defensive programming):
139 : mmaechler 2822 stop(gettextf("symmetry form '%s' is not yet implemented",
140 :     sym), domain = NA)
141 : maechler 2054 )
142 :     },
143 :     "pattern" = {
144 :     els <- scan(file, nmax = nz, quiet = TRUE,
145 : mmaechler 2369 what = list(i = integer(), j = integer()))
146 :     checkIJ(els)
147 : maechler 2054 switch(sym,
148 :     "general" = {
149 :     new("ngTMatrix", Dim = c(nr, nc),
150 :     i = els$i - 1L, j = els$j - 1L)
151 :     },
152 :     "symmetric" = {
153 :     new("nsTMatrix", uplo = "L", Dim = c(nr, nc),
154 :     i = els$i - 1L, j = els$j - 1L)
155 :     },
156 :     "skew-symmetric" = {
157 :     stop("symmetry form 'skew-symmetric' not yet implemented for reading")
158 :     ## FIXME: use dgT... but must expand the (i,j,x) slots!
159 :     new("ngTMatrix", uplo = "L", Dim = c(nr, nc),
160 :     i = els$i - 1L, j = els$j - 1L)
161 :    
162 :     },
163 :     "hermitian" = {
164 :     stop("symmetry form 'hermitian' not yet implemented for reading")
165 :     },
166 :     ## otherwise (not possible; just defensive programming):
167 : mmaechler 2822 stop(gettextf("symmetry form '%s' is not yet implemented",
168 :     sym), domain = NA)
169 : maechler 2054 )
170 :     },
171 :     "complex" = {
172 :     stop("element type 'complex' not yet implemented")
173 :     },
174 :     ## otherwise (not possible currently):
175 :     stop(gettextf("'%s()' is not yet implemented for element type '%s'",
176 : mmaechler 2822 "readMM", elt), domain = NA))
177 : maechler 2054 }
178 :     else
179 :     stop(gettextf("'%s()' is not yet implemented for representation '%s'",
180 : mmaechler 2822 "readMM", repr), domain = NA)
181 : bates 852 }

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