1 /* testExtract.c */
2
3 #include "../InpMtx.h"
4 #include "../../Drand.h"
5
6 /*--------------------------------------------------------------------*/
7 int
main(int argc,char * argv[])8 main ( int argc, char *argv[] )
9 /*
10 ------------------------------------------------------------------
11 generate a random matrix and test the submatrix extraction method.
12 the output is a matlab file to test correctness.
13
14 created -- 98oct15, cca
15 --------------------------------------------------------------------
16 */
17 {
18 int coordType, dataType, ii, msglvl, ncolA, ncol1, ncol2,
19 nitem, nrowA, nrow1, nrow2, rc, seed, symmetryflag ;
20 int *colids, *rowids, *temp ;
21 InpMtx *A ;
22 IV *cols1IV, *cols2IV, *rows1IV, *rows2IV ;
23 FILE *msgFile ;
24
25 if ( argc != 14 ) {
26 fprintf(stdout,
27 "\n\n %% usage : %s msglvl msgFile dataType symmetryflag coordType "
28 "\n %% nrowA ncolA nitem nrow1 nrow2 ncol1 ncol2 seed"
29 "\n %% msglvl -- message level"
30 "\n %% msgFile -- message file"
31 "\n %% dataType -- type of matrix entries"
32 "\n %% 1 -- real"
33 "\n %% 2 -- complex"
34 "\n %% symmetryflag -- type of matrix symmetry"
35 "\n %% 0 -- symmetric"
36 "\n %% 1 -- hermitian"
37 "\n %% 2 -- nonsymmetric"
38 "\n %% coordType -- coordinate Type"
39 "\n %% 1 -- by rows"
40 "\n %% 2 -- by columns"
41 "\n %% 3 -- by chevrons"
42 "\n %% nrowA -- number of rows in A"
43 "\n %% ncolA -- number of columns in A"
44 "\n %% nitem -- number of items to be loaded into A"
45 "\n %% nrow1 -- number of rows in the first block"
46 "\n %% nrow2 -- number of rows in the second block"
47 "\n %% ncol1 -- number of columns in the first block"
48 "\n %% ncol2 -- number of columns in the second block"
49 "\n %% seed -- random number seed"
50 "\n", argv[0]) ;
51 return(0) ;
52 }
53 msglvl = atoi(argv[1]) ;
54 if ( strcmp(argv[2], "stdout") == 0 ) {
55 msgFile = stdout ;
56 } else if ( (msgFile = fopen(argv[2], "a")) == NULL ) {
57 fprintf(stderr, "\n fatal error in %s"
58 "\n unable to open file %s\n",
59 argv[0], argv[2]) ;
60 return(-1) ;
61 }
62 dataType = atoi(argv[3]) ;
63 symmetryflag = atoi(argv[4]) ;
64 coordType = atoi(argv[5]) ;
65 nrowA = atoi(argv[6]) ;
66 ncolA = atoi(argv[7]) ;
67 nitem = atoi(argv[8]) ;
68 nrow1 = atoi(argv[9]) ;
69 nrow2 = atoi(argv[10]) ;
70 ncol1 = atoi(argv[11]) ;
71 ncol2 = atoi(argv[12]) ;
72 seed = atoi(argv[13]) ;
73 fprintf(msgFile,
74 "\n %% %s "
75 "\n %% msglvl -- %d"
76 "\n %% msgFile -- %s"
77 "\n %% dataType -- %d"
78 "\n %% symmetryflag -- %d"
79 "\n %% coordType -- %d"
80 "\n %% nrowA -- %d"
81 "\n %% ncolA -- %d"
82 "\n %% nitem -- %d"
83 "\n %% nrow1 -- %d"
84 "\n %% nrow2 -- %d"
85 "\n %% ncol1 -- %d"
86 "\n %% ncol2 -- %d"
87 "\n %% seed -- %d"
88 "\n",
89 argv[0], msglvl, argv[2], dataType, symmetryflag, coordType,
90 nrowA, ncolA, nitem, nrow1, nrow2, ncol1, ncol2, seed) ;
91 fflush(msgFile) ;
92 if ( dataType != SPOOLES_REAL && dataType != SPOOLES_COMPLEX ) {
93 fprintf(stderr, "\n invalid value %d for dataType\n", dataType) ;
94 exit(-1) ;
95 }
96 if ( symmetryflag != SPOOLES_SYMMETRIC
97 && symmetryflag != SPOOLES_HERMITIAN
98 && symmetryflag != SPOOLES_NONSYMMETRIC ) {
99 fprintf(stderr,
100 "\n invalid value %d for symmetryflag\n", symmetryflag) ;
101 exit(-1) ;
102 }
103 if ( symmetryflag == SPOOLES_HERMITIAN && dataType != SPOOLES_COMPLEX ){
104 fprintf(stderr,
105 "\n symmetryflag is hermitian, data type is not complex\n") ;
106 exit(-1) ;
107 }
108 if ( coordType != INPMTX_BY_ROWS
109 && coordType != INPMTX_BY_COLUMNS
110 && coordType != INPMTX_BY_CHEVRONS ) {
111 fprintf(stderr,
112 "\n invalid value %d for coordType\n", coordType) ;
113 exit(-1) ;
114 }
115 if ( nrowA <= 0 || ncolA <= 0 || nitem <= 0 ) {
116 fprintf(stderr,
117 "\n invalid value: nrow = %d, ncol = %d, nitem = %d",
118 nrowA, ncolA, nitem) ;
119 exit(-1) ;
120 }
121 fprintf(msgFile, "\n %% symmetryflag %d, nrowA %d, ncolA %d\n",
122 symmetryflag, nrowA, ncolA) ;
123 if ( ( symmetryflag == SPOOLES_SYMMETRIC
124 || symmetryflag == SPOOLES_HERMITIAN)
125 && nrowA != ncolA ) {
126 fprintf(stderr,
127 "\n symmetric matrix, nrowA %d, ncolA %d\n", nrowA, ncolA) ;
128 exit(-1) ;
129 }
130 if ( nrow1 < 0 || nrow2 < 0 || (nrow1 + nrow2 != nrowA) ) {
131 fprintf(stderr,
132 "\n invalid value: nrow = %d, nrow1 = %d, nrow2 = %d",
133 nrowA, nrow1, nrow2) ;
134 exit(-1) ;
135 }
136 if ( ncol1 < 0 || ncol2 < 0 || (ncol1 + ncol2 != ncolA) ) {
137 fprintf(stderr,
138 "\n invalid value: ncol = %d, ncol1 = %d, ncol2 = %d",
139 ncolA, ncol1, ncol2) ;
140 exit(-1) ;
141 }
142 /*
143 ----------------------------
144 initialize the matrix object
145 ----------------------------
146 */
147 A = InpMtx_new() ;
148 rc = InpMtx_randomMatrix(A, dataType, coordType, INPMTX_BY_VECTORS,
149 nrowA, ncolA, symmetryflag, 0, nitem, seed) ;
150 /*
151 -------------------------------------------
152 write the assembled matrix to a matlab file
153 -------------------------------------------
154 */
155 fprintf(msgFile, "\n A = zeros(%d,%d) ;", nrowA, ncolA) ;
156 InpMtx_writeForMatlab(A, "A", msgFile) ;
157 if ( symmetryflag == SPOOLES_SYMMETRIC ) {
158 fprintf(msgFile,
159 "\n A = diag(diag(A)) + triu(A,1) + transpose(triu(A,1)) ;") ;
160 } else if ( symmetryflag == SPOOLES_HERMITIAN ) {
161 fprintf(msgFile,
162 "\n A = diag(diag(A)) + triu(A,1) + ctranspose(triu(A,1)) ;") ;
163 }
164 /*
165 ----------------------------------
166 generate row and column id vectors
167 ----------------------------------
168 */
169 temp = IVinit(nrowA, -1) ;
170 IVramp(nrowA, temp, 0, 1) ;
171 IVshuffle(nrowA, temp, ++seed) ;
172 if ( nrow1 > 0 ) {
173 rows1IV = IV_new() ;
174 IV_init(rows1IV, nrow1, NULL) ;
175 IVqsortUp(nrow1, temp) ;
176 IVcopy(nrow1, IV_entries(rows1IV), temp) ;
177 } else {
178 rows1IV = NULL ;
179 }
180 if ( nrow2 > 0 ) {
181 rows2IV = IV_new() ;
182 IV_init(rows2IV, nrow2, NULL) ;
183 IVqsortUp(nrow2, temp + nrow1) ;
184 IVcopy(nrow2, IV_entries(rows2IV), temp + nrow1) ;
185 } else {
186 rows2IV = NULL ;
187 }
188 IVfree(temp) ;
189 if ( symmetryflag == SPOOLES_NONSYMMETRIC ) {
190 temp = IVinit(ncolA, -1) ;
191 IVramp(ncolA, temp, 0, 1) ;
192 IVshuffle(ncolA, temp, ++seed) ;
193 if ( ncol1 > 0 ) {
194 cols1IV = IV_new() ;
195 IV_init(cols1IV, ncol1, NULL) ;
196 IVqsortUp(ncol1, temp) ;
197 IVcopy(ncol1, IV_entries(cols1IV), temp) ;
198 } else {
199 cols1IV = NULL ;
200 }
201 if ( ncol2 > 0 ) {
202 cols2IV = IV_new() ;
203 IV_init(cols2IV, ncol2, NULL) ;
204 IVqsortUp(ncol2, temp + ncol1) ;
205 IVcopy(ncol2, IV_entries(cols2IV), temp + ncol1) ;
206 } else {
207 cols2IV = NULL ;
208 }
209 IVfree(temp) ;
210 } else {
211 if ( ncol1 > 0 ) {
212 cols1IV = IV_new() ;
213 IV_init(cols1IV, ncol1, NULL) ;
214 IV_copy(cols1IV, rows1IV) ;
215 } else {
216 cols1IV = NULL ;
217 }
218 if ( ncol2 > 0 ) {
219 cols2IV = IV_new() ;
220 IV_init(cols2IV, ncol2, NULL) ;
221 IV_copy(cols2IV, rows2IV) ;
222 } else {
223 cols2IV = NULL ;
224 }
225 }
226 if ( nrow1 > 0 ) {
227 fprintf(msgFile, "\n rows1 = zeros(%d,1) ;", nrow1) ;
228 IV_writeForMatlab(rows1IV, "rows1", msgFile) ;
229 }
230 if ( nrow2 > 0 ) {
231 fprintf(msgFile, "\n rows2 = zeros(%d,1) ;", nrow2) ;
232 IV_writeForMatlab(rows2IV, "rows2", msgFile) ;
233 }
234 if ( ncol1 > 0 ) {
235 fprintf(msgFile, "\n cols1 = zeros(%d,1) ;", ncol1) ;
236 IV_writeForMatlab(cols1IV, "cols1", msgFile) ;
237 }
238 if ( ncol2 > 0 ) {
239 fprintf(msgFile, "\n cols2 = zeros(%d,1) ;", ncol2) ;
240 IV_writeForMatlab(cols2IV, "cols2", msgFile) ;
241 }
242 /*
243 -----------------------
244 extract the submatrices
245 -----------------------
246 */
247 if ( nrow1 > 0 ) {
248 if ( ncol1 > 0 ) {
249 InpMtx *A11 = InpMtx_new() ;
250 rc = InpMtx_initFromSubmatrix(A11, A, rows1IV, cols1IV,
251 symmetryflag, msglvl, msgFile) ;
252 if ( rc != 1 ) {
253 fprintf(stderr, "\n error return %d for A11\n", rc) ;
254 exit(-1) ;
255 }
256 fprintf(msgFile, "\n A11 = zeros(%d,%d) ;", nrow1, ncol1) ;
257 InpMtx_writeForMatlab(A11, "A11", msgFile) ;
258 if ( symmetryflag == SPOOLES_SYMMETRIC ) {
259 fprintf(msgFile,
260 "\n A11 = diag(diag(A11)) + triu(A11,1) + transpose(triu(A11,1)) ;") ;
261 } else if ( symmetryflag == SPOOLES_HERMITIAN ) {
262 fprintf(msgFile,
263 "\n A11 = diag(diag(A11)) + triu(A11,1) + ctranspose(triu(A11,1)) ;") ;
264 }
265 fprintf(msgFile,
266 "\n err11 = max(max(abs(A11 - A(rows1,cols1)))) ;") ;
267 InpMtx_free(A11) ;
268 } else {
269 fprintf(msgFile, "\n err11 = 0 ;") ;
270 }
271 if ( ncol2 > 0 ) {
272 InpMtx *A12 = InpMtx_new() ;
273 rc = InpMtx_initFromSubmatrix(A12, A, rows1IV, cols2IV,
274 symmetryflag, msglvl, msgFile) ;
275 if ( rc != 1 ) {
276 fprintf(stderr, "\n error return %d for A12\n", rc) ;
277 exit(-1) ;
278 }
279 fprintf(msgFile, "\n A12 = zeros(%d,%d) ;", nrow1, ncol2) ;
280 InpMtx_writeForMatlab(A12, "A12", msgFile) ;
281 fprintf(msgFile,
282 "\n err12 = max(max(abs(A12 - A(rows1,cols2)))) ;") ;
283 InpMtx_free(A12) ;
284 } else {
285 fprintf(msgFile, "\n err12 = 0 ;") ;
286 }
287 } else {
288 fprintf(msgFile, "\n err11 = 0 ;") ;
289 fprintf(msgFile, "\n err12 = 0 ;") ;
290 }
291 if ( nrow2 > 0 ) {
292 if ( ncol1 > 0 ) {
293 InpMtx *A21 = InpMtx_new() ;
294 rc = InpMtx_initFromSubmatrix(A21, A, rows2IV, cols1IV,
295 symmetryflag, msglvl, msgFile) ;
296 if ( rc != 1 ) {
297 fprintf(stderr, "\n error return %d for A21\n", rc) ;
298 exit(-1) ;
299 }
300 fprintf(msgFile, "\n A21 = zeros(%d,%d) ;", nrow2, ncol1) ;
301 InpMtx_writeForMatlab(A21, "A21", msgFile) ;
302 fprintf(msgFile,
303 "\n err21 = max(max(abs(A21 - A(rows2,cols1)))) ;") ;
304 InpMtx_free(A21) ;
305 } else {
306 fprintf(msgFile, "\n err21 = 0 ;") ;
307 }
308 if ( ncol2 > 0 ) {
309 InpMtx *A22 = InpMtx_new() ;
310 rc = InpMtx_initFromSubmatrix(A22, A, rows2IV, cols2IV,
311 symmetryflag, msglvl, msgFile) ;
312 if ( rc != 1 ) {
313 fprintf(stderr, "\n error return %d for A22\n", rc) ;
314 exit(-1) ;
315 }
316 fprintf(msgFile, "\n A22 = zeros(%d,%d) ;", nrow2, ncol2) ;
317 InpMtx_writeForMatlab(A22, "A22", msgFile) ;
318 if ( symmetryflag == SPOOLES_SYMMETRIC ) {
319 fprintf(msgFile,
320 "\n A22 = diag(diag(A22)) + triu(A22,1) + transpose(triu(A22,1)) ;") ;
321 } else if ( symmetryflag == SPOOLES_HERMITIAN ) {
322 fprintf(msgFile,
323 "\n A22 = diag(diag(A22)) + triu(A22,1) + ctranspose(triu(A22,1)) ;") ;
324 }
325 fprintf(msgFile,
326 "\n err22 = max(max(abs(A22 - A(rows2,cols2)))) ;") ;
327 InpMtx_free(A22) ;
328 } else {
329 fprintf(msgFile, "\n err22 = 0 ;") ;
330 }
331 } else {
332 fprintf(msgFile, "\n err21 = 0 ;") ;
333 fprintf(msgFile, "\n err22 = 0 ;") ;
334 }
335 fprintf(msgFile, "\n error = [ err11 err12 ; err21 err22 ] ") ;
336 /*
337 ------------------------
338 free the working storage
339 ------------------------
340 */
341 InpMtx_free(A) ;
342 if ( rows1IV != NULL ) { IV_free(rows1IV) ; }
343 if ( rows2IV != NULL ) { IV_free(rows2IV) ; }
344 if ( cols1IV != NULL ) { IV_free(cols1IV) ; }
345 if ( cols2IV != NULL ) { IV_free(cols2IV) ; }
346
347 fclose(msgFile) ;
348
349 return(1) ; }
350
351 /*--------------------------------------------------------------------*/
352