SCM

SCM Repository

[matrix] Annotation of /branches/Matrix-mer2/src/iohb.c
ViewVC logotype

Annotation of /branches/Matrix-mer2/src/iohb.c

Parent Directory Parent Directory | Revision Log Revision Log


Revision 825 - (view) (download) (as text)
Original Path: pkg/src/iohb.c

1 : bates 825 /*
2 :     Fri Aug 15 16:29:47 EDT 1997
3 :    
4 :     Harwell-Boeing File I/O in C
5 :     V. 1.0
6 :    
7 :     National Institute of Standards and Technology, MD.
8 :     K.A. Remington
9 :    
10 :     ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
11 :     NOTICE
12 :    
13 :     Permission to use, copy, modify, and distribute this software and
14 :     its documentation for any purpose and without fee is hereby granted
15 :     provided that the above copyright notice appear in all copies and
16 :     that both the copyright notice and this permission notice appear in
17 :     supporting documentation.
18 :    
19 :     Neither the Author nor the Institution (National Institute of Standards
20 :     and Technology) make any representations about the suitability of this
21 :     software for any purpose. This software is provided "as is" without
22 :     expressed or implied warranty.
23 :     ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
24 :    
25 :     ---------------------
26 :     INTERFACE DESCRIPTION
27 :     ---------------------
28 :     ---------------
29 :     QUERY FUNCTIONS
30 :     ---------------
31 :    
32 :     FUNCTION:
33 :    
34 :     int readHB_info(const char *filename, int *M, int *N, int *nz,
35 :     char **Type, int *Nrhs)
36 :    
37 :     DESCRIPTION:
38 :    
39 :     The readHB_info function opens and reads the header information from
40 :     the specified Harwell-Boeing file, and reports back the number of rows
41 :     and columns in the stored matrix (M and N), the number of nonzeros in
42 :     the matrix (nz), the 3-character matrix type(Type), and the number of
43 :     right-hand-sides stored along with the matrix (Nrhs). This function
44 :     is designed to retrieve basic size information which can be used to
45 :     allocate arrays.
46 :    
47 :     FUNCTION:
48 :    
49 :     int readHB_header(FILE* in_file, char* Title, char* Key, char* Type,
50 :     int* Nrow, int* Ncol, int* Nnzero, int* Nrhs,
51 :     char* Ptrfmt, char* Indfmt, char* Valfmt, char* Rhsfmt,
52 :     int* Ptrcrd, int* Indcrd, int* Valcrd, int* Rhscrd,
53 :     char *Rhstype)
54 :    
55 :     DESCRIPTION:
56 :    
57 :     More detailed than the readHB_info function, readHB_header() reads from
58 :     the specified Harwell-Boeing file all of the header information.
59 :    
60 :    
61 :     ------------------------------
62 :     DOUBLE PRECISION I/O FUNCTIONS
63 :     ------------------------------
64 :     FUNCTION:
65 :    
66 :     int readHB_newmat_double(const char *filename, int *M, int *N, *int nz,
67 :     int **colptr, int **rowind, double**val)
68 :    
69 :     int readHB_mat_double(const char *filename, int *colptr, int *rowind,
70 :     double*val)
71 :    
72 :    
73 :     DESCRIPTION:
74 :    
75 :     This function opens and reads the specified file, interpreting its
76 :     contents as a sparse matrix stored in the Harwell/Boeing standard
77 :     format. (See readHB_aux_double to read auxillary vectors.)
78 :     -- Values are interpreted as double precision numbers. --
79 :    
80 :     The "mat" function uses _pre-allocated_ vectors to hold the index and
81 :     nonzero value information.
82 :    
83 :     The "newmat" function allocates vectors to hold the index and nonzero
84 :     value information, and returns pointers to these vectors along with
85 :     matrix dimension and number of nonzeros.
86 :    
87 :     FUNCTION:
88 :    
89 :     int readHB_aux_double(const char* filename, const char AuxType, double b[])
90 :    
91 :     int readHB_newaux_double(const char* filename, const char AuxType, double** b)
92 :    
93 :     DESCRIPTION:
94 :    
95 :     This function opens and reads from the specified file auxillary vector(s).
96 :     The char argument Auxtype determines which type of auxillary vector(s)
97 :     will be read (if present in the file).
98 :    
99 :     AuxType = 'F' right-hand-side
100 :     AuxType = 'G' initial estimate (Guess)
101 :     AuxType = 'X' eXact solution
102 :    
103 :     If Nrhs > 1, all of the Nrhs vectors of the given type are read and
104 :     stored in column-major order in the vector b.
105 :    
106 :     The "newaux" function allocates a vector to hold the values retrieved.
107 :     The "mat" function uses a _pre-allocated_ vector to hold the values.
108 :    
109 :     FUNCTION:
110 :    
111 :     int writeHB_mat_double(const char* filename, int M, int N,
112 :     int nz, const int colptr[], const int rowind[],
113 :     const double val[], int Nrhs, const double rhs[],
114 :     const double guess[], const double exact[],
115 :     const char* Title, const char* Key, const char* Type,
116 :     char* Ptrfmt, char* Indfmt, char* Valfmt, char* Rhsfmt,
117 :     const char* Rhstype)
118 :    
119 :     DESCRIPTION:
120 :    
121 :     The writeHB_mat_double function opens the named file and writes the specified
122 :     matrix and optional auxillary vector(s) to that file in Harwell-Boeing
123 :     format. The format arguments (Ptrfmt,Indfmt,Valfmt, and Rhsfmt) are
124 :     character strings specifying "Fortran-style" output formats -- as they
125 :     would appear in a Harwell-Boeing file. They are used to produce output
126 :     which is as close as possible to what would be produced by Fortran code,
127 :     but note that "D" and "P" edit descriptors are not supported.
128 :     If NULL, the following defaults will be used:
129 :     Ptrfmt = Indfmt = "(8I10)"
130 :     Valfmt = Rhsfmt = "(4E20.13)"
131 :    
132 :     -----------------------
133 :     CHARACTER I/O FUNCTIONS
134 :     -----------------------
135 :     FUNCTION:
136 :    
137 :     int readHB_mat_char(const char* filename, int colptr[], int rowind[],
138 :     char val[], char* Valfmt)
139 :     int readHB_newmat_char(const char* filename, int* M, int* N, int* nonzeros,
140 :     int** colptr, int** rowind, char** val, char** Valfmt)
141 :    
142 :     DESCRIPTION:
143 :    
144 :     This function opens and reads the specified file, interpreting its
145 :     contents as a sparse matrix stored in the Harwell/Boeing standard
146 :     format. (See readHB_aux_char to read auxillary vectors.)
147 :     -- Values are interpreted as char strings. --
148 :     (Used to translate exact values from the file into a new storage format.)
149 :    
150 :     The "mat" function uses _pre-allocated_ arrays to hold the index and
151 :     nonzero value information.
152 :    
153 :     The "newmat" function allocates char arrays to hold the index
154 :     and nonzero value information, and returns pointers to these arrays
155 :     along with matrix dimension and number of nonzeros.
156 :    
157 :     FUNCTION:
158 :    
159 :     int readHB_aux_char(const char* filename, const char AuxType, char b[])
160 :     int readHB_newaux_char(const char* filename, const char AuxType, char** b,
161 :     char** Rhsfmt)
162 :    
163 :     DESCRIPTION:
164 :    
165 :     This function opens and reads from the specified file auxillary vector(s).
166 :     The char argument Auxtype determines which type of auxillary vector(s)
167 :     will be read (if present in the file).
168 :    
169 :     AuxType = 'F' right-hand-side
170 :     AuxType = 'G' initial estimate (Guess)
171 :     AuxType = 'X' eXact solution
172 :    
173 :     If Nrhs > 1, all of the Nrhs vectors of the given type are read and
174 :     stored in column-major order in the vector b.
175 :    
176 :     The "newaux" function allocates a character array to hold the values
177 :     retrieved.
178 :     The "mat" function uses a _pre-allocated_ array to hold the values.
179 :    
180 :     FUNCTION:
181 :    
182 :     int writeHB_mat_char(const char* filename, int M, int N,
183 :     int nz, const int colptr[], const int rowind[],
184 :     const char val[], int Nrhs, const char rhs[],
185 :     const char guess[], const char exact[],
186 :     const char* Title, const char* Key, const char* Type,
187 :     char* Ptrfmt, char* Indfmt, char* Valfmt, char* Rhsfmt,
188 :     const char* Rhstype)
189 :    
190 :     DESCRIPTION:
191 :    
192 :     The writeHB_mat_char function opens the named file and writes the specified
193 :     matrix and optional auxillary vector(s) to that file in Harwell-Boeing
194 :     format. The format arguments (Ptrfmt,Indfmt,Valfmt, and Rhsfmt) are
195 :     character strings specifying "Fortran-style" output formats -- as they
196 :     would appear in a Harwell-Boeing file. Valfmt and Rhsfmt must accurately
197 :     represent the character representation of the values stored in val[]
198 :     and rhs[].
199 :    
200 :     If NULL, the following defaults will be used for the integer vectors:
201 :     Ptrfmt = Indfmt = "(8I10)"
202 :     Valfmt = Rhsfmt = "(4E20.13)"
203 :    
204 :    
205 :     */
206 :    
207 :     /*---------------------------------------------------------------------*/
208 :     /* If zero-based indexing is desired, _SP_base should be set to 0 */
209 :     /* This will cause indices read from H-B files to be decremented by 1 */
210 :     /* and indices written to H-B files to be incremented by 1 */
211 :     /* <<< Standard usage is _SP_base = 1 >>> */
212 :     /* <<< Changed to _SP_base = 0 DMB 2005-08-06 >>> */
213 :     #ifndef _SP_base
214 :     #define _SP_base 0
215 :     #endif
216 :     /*---------------------------------------------------------------------*/
217 :    
218 :     #include "iohb.h"
219 :     #include<stdio.h>
220 :     #include<stdlib.h>
221 :     #include<string.h>
222 :     #include<math.h>
223 :     #include<malloc.h>
224 :    
225 :     char* substr(const char* S, const int pos, const int len);
226 :     void upcase(char* S);
227 :     void IOHBTerminate(char* message);
228 :    
229 :     int readHB_info(const char* filename, int* M, int* N, int* nz, char** Type,
230 :     int* Nrhs)
231 :     {
232 :     /****************************************************************************/
233 :     /* The readHB_info function opens and reads the header information from */
234 :     /* the specified Harwell-Boeing file, and reports back the number of rows */
235 :     /* and columns in the stored matrix (M and N), the number of nonzeros in */
236 :     /* the matrix (nz), and the number of right-hand-sides stored along with */
237 :     /* the matrix (Nrhs). */
238 :     /* */
239 :     /* For a description of the Harwell Boeing standard, see: */
240 :     /* Duff, et al., ACM TOMS Vol.15, No.1, March 1989 */
241 :     /* */
242 :     /* ---------- */
243 :     /* **CAVEAT** */
244 :     /* ---------- */
245 :     /* ** If the input file does not adhere to the H/B format, the ** */
246 :     /* ** results will be unpredictable. ** */
247 :     /* */
248 :     /****************************************************************************/
249 :     FILE *in_file;
250 :     int Ptrcrd, Indcrd, Valcrd, Rhscrd;
251 :     int Nrow, Ncol, Nnzero;
252 :     char *mat_type;
253 :     char Title[73], Key[9], Rhstype[4];
254 :     char Ptrfmt[17], Indfmt[17], Valfmt[21], Rhsfmt[21];
255 :    
256 :     mat_type = (char *) malloc(4);
257 :     if ( mat_type == NULL ) IOHBTerminate("Insufficient memory for mat_typen");
258 :    
259 :     if ( (in_file = fopen( filename, "r")) == NULL ) {
260 :     fprintf(stderr,"Error: Cannot open file: %s\n",filename);
261 :     return 0;
262 :     }
263 :    
264 :     readHB_header(in_file, Title, Key, mat_type, &Nrow, &Ncol, &Nnzero, Nrhs,
265 :     Ptrfmt, Indfmt, Valfmt, Rhsfmt,
266 :     &Ptrcrd, &Indcrd, &Valcrd, &Rhscrd, Rhstype);
267 :     fclose(in_file);
268 :     *Type = mat_type;
269 :     *(*Type+3) = (char) NULL;
270 :     *M = Nrow;
271 :     *N = Ncol;
272 :     *nz = Nnzero;
273 :     if (Rhscrd == 0) {*Nrhs = 0;}
274 :    
275 :     /* In verbose mode, print some of the header information: */
276 :     /*
277 :     if (verbose == 1)
278 :     {
279 :     printf("Reading from Harwell-Boeing file %s (verbose on)...\n",filename);
280 :     printf(" Title: %s\n",Title);
281 :     printf(" Key: %s\n",Key);
282 :     printf(" The stored matrix is %i by %i with %i nonzeros.\n",
283 :     *M, *N, *nz );
284 :     printf(" %i right-hand--side(s) stored.\n",*Nrhs);
285 :     }
286 :     */
287 :    
288 :     return 1;
289 :    
290 :     }
291 :    
292 :    
293 :    
294 :     int readHB_header(FILE* in_file, char* Title, char* Key, char* Type,
295 :     int* Nrow, int* Ncol, int* Nnzero, int* Nrhs,
296 :     char* Ptrfmt, char* Indfmt, char* Valfmt, char* Rhsfmt,
297 :     int* Ptrcrd, int* Indcrd, int* Valcrd, int* Rhscrd,
298 :     char *Rhstype)
299 :     {
300 :     /*************************************************************************/
301 :     /* Read header information from the named H/B file... */
302 :     /*************************************************************************/
303 :     int Totcrd,Neltvl,Nrhsix;
304 :     char line[BUFSIZ];
305 :    
306 :     /* First line: */
307 :     fgets(line, BUFSIZ, in_file);
308 :     if ( sscanf(line,"%*s") < 0 )
309 :     IOHBTerminate("iohb.c: Null (or blank) first line of HB file.\n");
310 :     (void) sscanf(line, "%72c%8[^\n]", Title, Key);
311 :     *(Key+8) = (char) NULL;
312 :     *(Title+72) = (char) NULL;
313 :    
314 :     /* Second line: */
315 :     fgets(line, BUFSIZ, in_file);
316 :     if ( sscanf(line,"%*s") < 0 )
317 :     IOHBTerminate("iohb.c: Null (or blank) second line of HB file.\n");
318 :     if ( sscanf(line,"%i",&Totcrd) != 1) Totcrd = 0;
319 :     if ( sscanf(line,"%*i%i",Ptrcrd) != 1) *Ptrcrd = 0;
320 :     if ( sscanf(line,"%*i%*i%i",Indcrd) != 1) *Indcrd = 0;
321 :     if ( sscanf(line,"%*i%*i%*i%i",Valcrd) != 1) *Valcrd = 0;
322 :     if ( sscanf(line,"%*i%*i%*i%*i%i",Rhscrd) != 1) *Rhscrd = 0;
323 :    
324 :     /* Third line: */
325 :     fgets(line, BUFSIZ, in_file);
326 :     if ( sscanf(line,"%*s") < 0 )
327 :     IOHBTerminate("iohb.c: Null (or blank) third line of HB file.\n");
328 :     if ( sscanf(line, "%3c", Type) != 1)
329 :     IOHBTerminate("iohb.c: Invalid Type info, line 3 of Harwell-Boeing file.\n");
330 :     upcase(Type);
331 :     if ( sscanf(line,"%*3c%i",Nrow) != 1) *Nrow = 0 ;
332 :     if ( sscanf(line,"%*3c%*i%i",Ncol) != 1) *Ncol = 0 ;
333 :     if ( sscanf(line,"%*3c%*i%*i%i",Nnzero) != 1) *Nnzero = 0 ;
334 :     if ( sscanf(line,"%*3c%*i%*i%*i%i",&Neltvl) != 1) Neltvl = 0 ;
335 :    
336 :     /* Fourth line: */
337 :     fgets(line, BUFSIZ, in_file);
338 :     if ( sscanf(line,"%*s") < 0 )
339 :     IOHBTerminate("iohb.c: Null (or blank) fourth line of HB file.\n");
340 :     if ( sscanf(line, "%16c",Ptrfmt) != 1)
341 :     IOHBTerminate("iohb.c: Invalid format info, line 4 of Harwell-Boeing file.\n");
342 :     if ( sscanf(line, "%*16c%16c",Indfmt) != 1)
343 :     IOHBTerminate("iohb.c: Invalid format info, line 4 of Harwell-Boeing file.\n");
344 :     if ( sscanf(line, "%*16c%*16c%20c",Valfmt) != 1)
345 :     IOHBTerminate("iohb.c: Invalid format info, line 4 of Harwell-Boeing file.\n");
346 :     sscanf(line, "%*16c%*16c%*20c%20c",Rhsfmt);
347 :     *(Ptrfmt+16) = (char) NULL;
348 :     *(Indfmt+16) = (char) NULL;
349 :     *(Valfmt+20) = (char) NULL;
350 :     *(Rhsfmt+20) = (char) NULL;
351 :    
352 :     /* (Optional) Fifth line: */
353 :     if (*Rhscrd != 0 )
354 :     {
355 :     fgets(line, BUFSIZ, in_file);
356 :     if ( sscanf(line,"%*s") < 0 )
357 :     IOHBTerminate("iohb.c: Null (or blank) fifth line of HB file.\n");
358 :     if ( sscanf(line, "%3c", Rhstype) != 1)
359 :     IOHBTerminate("iohb.c: Invalid RHS type information, line 5 of Harwell-Boeing file.\n");
360 :     if ( sscanf(line, "%*3c%i", Nrhs) != 1) *Nrhs = 0;
361 :     if ( sscanf(line, "%*3c%*i%i", &Nrhsix) != 1) Nrhsix = 0;
362 :     }
363 :     return 1;
364 :     }
365 :    
366 :    
367 :     int readHB_mat_double(const char* filename, int colptr[], int rowind[],
368 :     double val[])
369 :     {
370 :     /****************************************************************************/
371 :     /* This function opens and reads the specified file, interpreting its */
372 :     /* contents as a sparse matrix stored in the Harwell/Boeing standard */
373 :     /* format and creating compressed column storage scheme vectors to hold */
374 :     /* the index and nonzero value information. */
375 :     /* */
376 :     /* ---------- */
377 :     /* **CAVEAT** */
378 :     /* ---------- */
379 :     /* Parsing real formats from Fortran is tricky, and this file reader */
380 :     /* does not claim to be foolproof. It has been tested for cases when */
381 :     /* the real values are printed consistently and evenly spaced on each */
382 :     /* line, with Fixed (F), and Exponential (E or D) formats. */
383 :     /* */
384 :     /* ** If the input file does not adhere to the H/B format, the ** */
385 :     /* ** results will be unpredictable. ** */
386 :     /* */
387 :     /****************************************************************************/
388 :     FILE *in_file;
389 :     int i,j,ind,col,offset,count,last,Nrhs;
390 :     int Ptrcrd, Indcrd, Valcrd, Rhscrd;
391 :     int Nrow, Ncol, Nnzero, Nentries;
392 :     int Ptrperline, Ptrwidth, Indperline, Indwidth;
393 :     int Valperline, Valwidth, Valprec;
394 :     int Valflag; /* Indicates 'E','D', or 'F' float format */
395 :     char* ThisElement;
396 :     char Title[73], Key[8], Type[4], Rhstype[4];
397 :     char Ptrfmt[17], Indfmt[17], Valfmt[21], Rhsfmt[21];
398 :     char line[BUFSIZ];
399 :    
400 :     if ( (in_file = fopen( filename, "r")) == NULL ) {
401 :     fprintf(stderr,"Error: Cannot open file: %s\n",filename);
402 :     return 0;
403 :     }
404 :    
405 :     readHB_header(in_file, Title, Key, Type, &Nrow, &Ncol, &Nnzero, &Nrhs,
406 :     Ptrfmt, Indfmt, Valfmt, Rhsfmt,
407 :     &Ptrcrd, &Indcrd, &Valcrd, &Rhscrd, Rhstype);
408 :    
409 :     /* Parse the array input formats from Line 3 of HB file */
410 :     ParseIfmt(Ptrfmt,&Ptrperline,&Ptrwidth);
411 :     ParseIfmt(Indfmt,&Indperline,&Indwidth);
412 :     if ( Type[0] != 'P' ) { /* Skip if pattern only */
413 :     ParseRfmt(Valfmt,&Valperline,&Valwidth,&Valprec,&Valflag);
414 :     }
415 :    
416 :     /* Read column pointer array: */
417 :    
418 :     offset = 1-_SP_base; /* if base 0 storage is declared (via macro definition), */
419 :     /* then storage entries are offset by 1 */
420 :    
421 :     ThisElement = (char *) malloc(Ptrwidth+1);
422 :     if ( ThisElement == NULL ) IOHBTerminate("Insufficient memory for ThisElement.");
423 :     *(ThisElement+Ptrwidth) = (char) NULL;
424 :     count=0;
425 :     for (i=0;i<Ptrcrd;i++)
426 :     {
427 :     fgets(line, BUFSIZ, in_file);
428 :     if ( sscanf(line,"%*s") < 0 )
429 :     IOHBTerminate("iohb.c: Null (or blank) line in pointer data region of HB file.\n");
430 :     col = 0;
431 :     for (ind = 0;ind<Ptrperline;ind++)
432 :     {
433 :     if (count > Ncol) break;
434 :     strncpy(ThisElement,line+col,Ptrwidth);
435 :     /* ThisElement = substr(line,col,Ptrwidth); */
436 :     colptr[count] = atoi(ThisElement)-offset;
437 :     count++; col += Ptrwidth;
438 :     }
439 :     }
440 :     free(ThisElement);
441 :    
442 :     /* Read row index array: */
443 :    
444 :     ThisElement = (char *) malloc(Indwidth+1);
445 :     if ( ThisElement == NULL ) IOHBTerminate("Insufficient memory for ThisElement.");
446 :     *(ThisElement+Indwidth) = (char) NULL;
447 :     count = 0;
448 :     for (i=0;i<Indcrd;i++)
449 :     {
450 :     fgets(line, BUFSIZ, in_file);
451 :     if ( sscanf(line,"%*s") < 0 )
452 :     IOHBTerminate("iohb.c: Null (or blank) line in index data region of HB file.\n");
453 :     col = 0;
454 :     for (ind = 0;ind<Indperline;ind++)
455 :     {
456 :     if (count == Nnzero) break;
457 :     strncpy(ThisElement,line+col,Indwidth);
458 :     /* ThisElement = substr(line,col,Indwidth); */
459 :     rowind[count] = atoi(ThisElement)-offset;
460 :     count++; col += Indwidth;
461 :     }
462 :     }
463 :     free(ThisElement);
464 :    
465 :     /* Read array of values: */
466 :    
467 :     if ( Type[0] != 'P' ) { /* Skip if pattern only */
468 :    
469 :     if ( Type[0] == 'C' ) Nentries = 2*Nnzero;
470 :     else Nentries = Nnzero;
471 :    
472 :     ThisElement = (char *) malloc(Valwidth+1);
473 :     if ( ThisElement == NULL ) IOHBTerminate("Insufficient memory for ThisElement.");
474 :     *(ThisElement+Valwidth) = (char) NULL;
475 :     count = 0;
476 :     for (i=0;i<Valcrd;i++)
477 :     {
478 :     fgets(line, BUFSIZ, in_file);
479 :     if ( sscanf(line,"%*s") < 0 )
480 :     IOHBTerminate("iohb.c: Null (or blank) line in value data region of HB file.\n");
481 :     if (Valflag == 'D') {
482 :     while( strchr(line,'D') ) *strchr(line,'D') = 'E';
483 :     /* *strchr(Valfmt,'D') = 'E'; */
484 :     }
485 :     col = 0;
486 :     for (ind = 0;ind<Valperline;ind++)
487 :     {
488 :     if (count == Nentries) break;
489 :     strncpy(ThisElement,line+col,Valwidth);
490 :     /*ThisElement = substr(line,col,Valwidth);*/
491 :     if ( Valflag != 'F' && strchr(ThisElement,'E') == NULL ) {
492 :     /* insert a char prefix for exp */
493 :     last = strlen(ThisElement);
494 :     for (j=last+1;j>=0;j--) {
495 :     ThisElement[j] = ThisElement[j-1];
496 :     if ( ThisElement[j] == '+' || ThisElement[j] == '-' ) {
497 :     ThisElement[j-1] = Valflag;
498 :     break;
499 :     }
500 :     }
501 :     }
502 :     val[count] = atof(ThisElement);
503 :     count++; col += Valwidth;
504 :     }
505 :     }
506 :     free(ThisElement);
507 :     }
508 :    
509 :     fclose(in_file);
510 :     return 1;
511 :     }
512 :    
513 :     int readHB_newmat_double(const char* filename, int* M, int* N, int* nonzeros,
514 :     int** colptr, int** rowind, double** val)
515 :     {
516 :     int Nrhs;
517 :     char *Type;
518 :    
519 :     readHB_info(filename, M, N, nonzeros, &Type, &Nrhs);
520 :    
521 :     *colptr = (int *)malloc((*N+1)*sizeof(int));
522 :     if ( *colptr == NULL ) IOHBTerminate("Insufficient memory for colptr.\n");
523 :     *rowind = (int *)malloc(*nonzeros*sizeof(int));
524 :     if ( *rowind == NULL ) IOHBTerminate("Insufficient memory for rowind.\n");
525 :     if ( Type[0] == 'C' ) {
526 :     /*
527 :     fprintf(stderr, "Warning: Reading complex data from HB file %s.\n",filename);
528 :     fprintf(stderr, " Real and imaginary parts will be interlaced in val[].\n");
529 :     */
530 :     /* Malloc enough space for real AND imaginary parts of val[] */
531 :     *val = (double *)malloc(*nonzeros*sizeof(double)*2);
532 :     if ( *val == NULL ) IOHBTerminate("Insufficient memory for val.\n");
533 :     } else {
534 :     if ( Type[0] != 'P' ) {
535 :     /* Malloc enough space for real array val[] */
536 :     *val = (double *)malloc(*nonzeros*sizeof(double));
537 :     if ( *val == NULL ) IOHBTerminate("Insufficient memory for val.\n");
538 :     }
539 :     } /* No val[] space needed if pattern only */
540 :     return readHB_mat_double(filename, *colptr, *rowind, *val);
541 :    
542 :     }
543 :    
544 :     int readHB_aux_double(const char* filename, const char AuxType, double b[])
545 :     {
546 :     /****************************************************************************/
547 :     /* This function opens and reads the specified file, placing auxillary */
548 :     /* vector(s) of the given type (if available) in b. */
549 :     /* Return value is the number of vectors successfully read. */
550 :     /* */
551 :     /* AuxType = 'F' full right-hand-side vector(s) */
552 :     /* AuxType = 'G' initial Guess vector(s) */
553 :     /* AuxType = 'X' eXact solution vector(s) */
554 :     /* */
555 :     /* ---------- */
556 :     /* **CAVEAT** */
557 :     /* ---------- */
558 :     /* Parsing real formats from Fortran is tricky, and this file reader */
559 :     /* does not claim to be foolproof. It has been tested for cases when */
560 :     /* the real values are printed consistently and evenly spaced on each */
561 :     /* line, with Fixed (F), and Exponential (E or D) formats. */
562 :     /* */
563 :     /* ** If the input file does not adhere to the H/B format, the ** */
564 :     /* ** results will be unpredictable. ** */
565 :     /* */
566 :     /****************************************************************************/
567 :     FILE *in_file;
568 :     int i,j,n,maxcol,start,stride,col,last,linel;
569 :     int Ptrcrd, Indcrd, Valcrd, Rhscrd;
570 :     int Nrow, Ncol, Nnzero, Nentries;
571 :     int Nrhs, nvecs, rhsi;
572 :     int Rhsperline, Rhswidth, Rhsprec;
573 :     int Rhsflag;
574 :     char *ThisElement;
575 :     char Title[73], Key[9], Type[4], Rhstype[4];
576 :     char Ptrfmt[17], Indfmt[17], Valfmt[21], Rhsfmt[21];
577 :     char line[BUFSIZ];
578 :    
579 :     if ((in_file = fopen( filename, "r")) == NULL) {
580 :     fprintf(stderr,"Error: Cannot open file: %s\n",filename);
581 :     return 0;
582 :     }
583 :    
584 :     readHB_header(in_file, Title, Key, Type, &Nrow, &Ncol, &Nnzero, &Nrhs,
585 :     Ptrfmt, Indfmt, Valfmt, Rhsfmt,
586 :     &Ptrcrd, &Indcrd, &Valcrd, &Rhscrd, Rhstype);
587 :    
588 :     if (Nrhs <= 0)
589 :     {
590 :     fprintf(stderr, "Warn: Attempt to read auxillary vector(s) when none are present.\n");
591 :     return 0;
592 :     }
593 :     if (Rhstype[0] != 'F' )
594 :     {
595 :     fprintf(stderr,"Warn: Attempt to read auxillary vector(s) which are not stored in Full form.\n");
596 :     fprintf(stderr," Rhs must be specified as full. \n");
597 :     return 0;
598 :     }
599 :    
600 :     /* If reading complex data, allow for interleaved real and imaginary values. */
601 :     if ( Type[0] == 'C' ) {
602 :     Nentries = 2*Nrow;
603 :     } else {
604 :     Nentries = Nrow;
605 :     }
606 :    
607 :     nvecs = 1;
608 :    
609 :     if ( Rhstype[1] == 'G' ) nvecs++;
610 :     if ( Rhstype[2] == 'X' ) nvecs++;
611 :    
612 :     if ( AuxType == 'G' && Rhstype[1] != 'G' ) {
613 :     fprintf(stderr, "Warn: Attempt to read auxillary Guess vector(s) when none are present.\n");
614 :     return 0;
615 :     }
616 :     if ( AuxType == 'X' && Rhstype[2] != 'X' ) {
617 :     fprintf(stderr, "Warn: Attempt to read auxillary eXact solution vector(s) when none are present.\n");
618 :     return 0;
619 :     }
620 :    
621 :     ParseRfmt(Rhsfmt, &Rhsperline, &Rhswidth, &Rhsprec,&Rhsflag);
622 :     maxcol = Rhsperline*Rhswidth;
623 :    
624 :     /* Lines to skip before starting to read RHS values... */
625 :     n = Ptrcrd + Indcrd + Valcrd;
626 :    
627 :     for (i = 0; i < n; i++)
628 :     fgets(line, BUFSIZ, in_file);
629 :    
630 :     /* start - number of initial aux vector entries to skip */
631 :     /* to reach first vector requested */
632 :     /* stride - number of aux vector entries to skip between */
633 :     /* requested vectors */
634 :     if ( AuxType == 'F' ) start = 0;
635 :     else if ( AuxType == 'G' ) start = Nentries;
636 :     else start = (nvecs-1)*Nentries;
637 :     stride = (nvecs-1)*Nentries;
638 :    
639 :     fgets(line, BUFSIZ, in_file);
640 :     linel= strchr(line,'\n')-line;
641 :     col = 0;
642 :     /* Skip to initial offset */
643 :    
644 :     for (i=0;i<start;i++) {
645 :     if ( col >= ( maxcol<linel?maxcol:linel ) ) {
646 :     fgets(line, BUFSIZ, in_file);
647 :     linel= strchr(line,'\n')-line;
648 :     col = 0;
649 :     }
650 :     col += Rhswidth;
651 :     }
652 :     if (Rhsflag == 'D') {
653 :     while( strchr(line,'D') ) *strchr(line,'D') = 'E';
654 :     }
655 :    
656 :     /* Read a vector of desired type, then skip to next */
657 :     /* repeating to fill Nrhs vectors */
658 :    
659 :     ThisElement = (char *) malloc(Rhswidth+1);
660 :     if ( ThisElement == NULL ) IOHBTerminate("Insufficient memory for ThisElement.");
661 :     *(ThisElement+Rhswidth) = (char) NULL;
662 :     for (rhsi=0;rhsi<Nrhs;rhsi++) {
663 :    
664 :     for (i=0;i<Nentries;i++) {
665 :     if ( col >= ( maxcol<linel?maxcol:linel ) ) {
666 :     fgets(line, BUFSIZ, in_file);
667 :     linel= strchr(line,'\n')-line;
668 :     if (Rhsflag == 'D') {
669 :     while( strchr(line,'D') ) *strchr(line,'D') = 'E';
670 :     }
671 :     col = 0;
672 :     }
673 :     strncpy(ThisElement,line+col,Rhswidth);
674 :     /*ThisElement = substr(line, col, Rhswidth);*/
675 :     if ( Rhsflag != 'F' && strchr(ThisElement,'E') == NULL ) {
676 :     /* insert a char prefix for exp */
677 :     last = strlen(ThisElement);
678 :     for (j=last+1;j>=0;j--) {
679 :     ThisElement[j] = ThisElement[j-1];
680 :     if ( ThisElement[j] == '+' || ThisElement[j] == '-' ) {
681 :     ThisElement[j-1] = Rhsflag;
682 :     break;
683 :     }
684 :     }
685 :     }
686 :     b[i] = atof(ThisElement);
687 :     col += Rhswidth;
688 :     }
689 :    
690 :     /* Skip any interleaved Guess/eXact vectors */
691 :    
692 :     for (i=0;i<stride;i++) {
693 :     if ( col >= ( maxcol<linel?maxcol:linel ) ) {
694 :     fgets(line, BUFSIZ, in_file);
695 :     linel= strchr(line,'\n')-line;
696 :     col = 0;
697 :     }
698 :     col += Rhswidth;
699 :     }
700 :    
701 :     }
702 :     free(ThisElement);
703 :    
704 :    
705 :     fclose(in_file);
706 :     return Nrhs;
707 :     }
708 :    
709 :     int readHB_newaux_double(const char* filename, const char AuxType, double** b)
710 :     {
711 :     int Nrhs,M,N,nonzeros;
712 :     char *Type;
713 :    
714 :     readHB_info(filename, &M, &N, &nonzeros, &Type, &Nrhs);
715 :     if ( Nrhs <= 0 ) {
716 :     fprintf(stderr,"Warn: Requested read of aux vector(s) when none are present.\n");
717 :     return 0;
718 :     } else {
719 :     if ( Type[0] == 'C' ) {
720 :     fprintf(stderr, "Warning: Reading complex aux vector(s) from HB file %s.",filename);
721 :     fprintf(stderr, " Real and imaginary parts will be interlaced in b[].");
722 :     *b = (double *)malloc(M*Nrhs*sizeof(double)*2);
723 :     if ( *b == NULL ) IOHBTerminate("Insufficient memory for rhs.\n");
724 :     return readHB_aux_double(filename, AuxType, *b);
725 :     } else {
726 :     *b = (double *)malloc(M*Nrhs*sizeof(double));
727 :     if ( *b == NULL ) IOHBTerminate("Insufficient memory for rhs.\n");
728 :     return readHB_aux_double(filename, AuxType, *b);
729 :     }
730 :     }
731 :     }
732 :    
733 :     int writeHB_mat_double(const char* filename, int M, int N,
734 :     int nz, const int colptr[], const int rowind[],
735 :     const double val[], int Nrhs, const double rhs[],
736 :     const double guess[], const double exact[],
737 :     const char* Title, const char* Key, const char* Type,
738 :     char* Ptrfmt, char* Indfmt, char* Valfmt, char* Rhsfmt,
739 :     const char* Rhstype)
740 :     {
741 :     /****************************************************************************/
742 :     /* The writeHB function opens the named file and writes the specified */
743 :     /* matrix and optional right-hand-side(s) to that file in Harwell-Boeing */
744 :     /* format. */
745 :     /* */
746 :     /* For a description of the Harwell Boeing standard, see: */
747 :     /* Duff, et al., ACM TOMS Vol.15, No.1, March 1989 */
748 :     /* */
749 :     /****************************************************************************/
750 :     FILE *out_file;
751 :     int i,j,entry,offset,acount,linemod;
752 :     int totcrd, ptrcrd, indcrd, valcrd, rhscrd;
753 :     int nvalentries, nrhsentries;
754 :     int Ptrperline, Ptrwidth, Indperline, Indwidth;
755 :     int Rhsperline, Rhswidth, Rhsprec;
756 :     int Rhsflag;
757 :     int Valperline, Valwidth, Valprec;
758 :     int Valflag; /* Indicates 'E','D', or 'F' float format */
759 :     char pformat[16],iformat[16],vformat[19],rformat[19];
760 :    
761 :     if ( Type[0] == 'C' ) {
762 :     nvalentries = 2*nz;
763 :     nrhsentries = 2*M;
764 :     } else {
765 :     nvalentries = nz;
766 :     nrhsentries = M;
767 :     }
768 :    
769 :     if ( filename != NULL ) {
770 :     if ( (out_file = fopen( filename, "w")) == NULL ) {
771 :     fprintf(stderr,"Error: Cannot open file: %s\n",filename);
772 :     return 0;
773 :     }
774 :     } else out_file = stdout;
775 :    
776 :     if ( Ptrfmt == NULL ) Ptrfmt = "(8I10)";
777 :     ParseIfmt(Ptrfmt,&Ptrperline,&Ptrwidth);
778 :     sprintf(pformat,"%%%dd",Ptrwidth);
779 :     ptrcrd = (N+1)/Ptrperline;
780 :     if ( (N+1)%Ptrperline != 0) ptrcrd++;
781 :    
782 :     if ( Indfmt == NULL ) Indfmt = Ptrfmt;
783 :     ParseIfmt(Indfmt,&Indperline,&Indwidth);
784 :     sprintf(iformat,"%%%dd",Indwidth);
785 :     indcrd = nz/Indperline;
786 :     if ( nz%Indperline != 0) indcrd++;
787 :    
788 :     if ( Type[0] != 'P' ) { /* Skip if pattern only */
789 :     if ( Valfmt == NULL ) Valfmt = "(4E20.13)";
790 :     ParseRfmt(Valfmt,&Valperline,&Valwidth,&Valprec,&Valflag);
791 :     if (Valflag == 'D') *strchr(Valfmt,'D') = 'E';
792 :     if (Valflag == 'F')
793 :     sprintf(vformat,"%% %d.%df",Valwidth,Valprec);
794 :     else
795 :     sprintf(vformat,"%% %d.%dE",Valwidth,Valprec);
796 :     valcrd = nvalentries/Valperline;
797 :     if ( nvalentries%Valperline != 0) valcrd++;
798 :     } else valcrd = 0;
799 :    
800 :     if ( Nrhs > 0 ) {
801 :     if ( Rhsfmt == NULL ) Rhsfmt = Valfmt;
802 :     ParseRfmt(Rhsfmt,&Rhsperline,&Rhswidth,&Rhsprec, &Rhsflag);
803 :     if (Rhsflag == 'F')
804 :     sprintf(rformat,"%% %d.%df",Rhswidth,Rhsprec);
805 :     else
806 :     sprintf(rformat,"%% %d.%dE",Rhswidth,Rhsprec);
807 :     if (Rhsflag == 'D') *strchr(Rhsfmt,'D') = 'E';
808 :     rhscrd = nrhsentries/Rhsperline;
809 :     if ( nrhsentries%Rhsperline != 0) rhscrd++;
810 :     if ( Rhstype[1] == 'G' ) rhscrd+=rhscrd;
811 :     if ( Rhstype[2] == 'X' ) rhscrd+=rhscrd;
812 :     rhscrd*=Nrhs;
813 :     } else rhscrd = 0;
814 :    
815 :     totcrd = 4+ptrcrd+indcrd+valcrd+rhscrd;
816 :    
817 :    
818 :     /* Print header information: */
819 :    
820 :     fprintf(out_file,"%-72s%-8s\n%14d%14d%14d%14d%14d\n",Title, Key, totcrd,
821 :     ptrcrd, indcrd, valcrd, rhscrd);
822 :     fprintf(out_file,"%3s%11s%14d%14d%14d\n",Type," ", M, N, nz);
823 :     fprintf(out_file,"%-16s%-16s%-20s", Ptrfmt, Indfmt, Valfmt);
824 :     if ( Nrhs != 0 ) {
825 :     /* Print Rhsfmt on fourth line and */
826 :     /* optional fifth header line for auxillary vector information: */
827 :     fprintf(out_file,"%-20s\n%-14s%d\n",Rhsfmt,Rhstype,Nrhs);
828 :     } else fprintf(out_file,"\n");
829 :    
830 :     offset = 1-_SP_base; /* if base 0 storage is declared (via macro definition), */
831 :     /* then storage entries are offset by 1 */
832 :    
833 :     /* Print column pointers: */
834 :     for (i=0;i<N+1;i++)
835 :     {
836 :     entry = colptr[i]+offset;
837 :     fprintf(out_file,pformat,entry);
838 :     if ( (i+1)%Ptrperline == 0 ) fprintf(out_file,"\n");
839 :     }
840 :    
841 :     if ( (N+1) % Ptrperline != 0 ) fprintf(out_file,"\n");
842 :    
843 :     /* Print row indices: */
844 :     for (i=0;i<nz;i++)
845 :     {
846 :     entry = rowind[i]+offset;
847 :     fprintf(out_file,iformat,entry);
848 :     if ( (i+1)%Indperline == 0 ) fprintf(out_file,"\n");
849 :     }
850 :    
851 :     if ( nz % Indperline != 0 ) fprintf(out_file,"\n");
852 :    
853 :     /* Print values: */
854 :    
855 :     if ( Type[0] != 'P' ) { /* Skip if pattern only */
856 :    
857 :     for (i=0;i<nvalentries;i++)
858 :     {
859 :     fprintf(out_file,vformat,val[i]);
860 :     if ( (i+1)%Valperline == 0 ) fprintf(out_file,"\n");
861 :     }
862 :    
863 :     if ( nvalentries % Valperline != 0 ) fprintf(out_file,"\n");
864 :    
865 :     /* If available, print right hand sides,
866 :     guess vectors and exact solution vectors: */
867 :     acount = 1;
868 :     linemod = 0;
869 :     if ( Nrhs > 0 ) {
870 :     for (i=0;i<Nrhs;i++)
871 :     {
872 :     for ( j=0;j<nrhsentries;j++ ) {
873 :     fprintf(out_file,rformat,rhs[j]);
874 :     if ( acount++%Rhsperline == linemod ) fprintf(out_file,"\n");
875 :     }
876 :     if ( acount%Rhsperline != linemod ) {
877 :     fprintf(out_file,"\n");
878 :     linemod = (acount-1)%Rhsperline;
879 :     }
880 :     rhs += nrhsentries;
881 :     if ( Rhstype[1] == 'G' ) {
882 :     for ( j=0;j<nrhsentries;j++ ) {
883 :     fprintf(out_file,rformat,guess[j]);
884 :     if ( acount++%Rhsperline == linemod ) fprintf(out_file,"\n");
885 :     }
886 :     if ( acount%Rhsperline != linemod ) {
887 :     fprintf(out_file,"\n");
888 :     linemod = (acount-1)%Rhsperline;
889 :     }
890 :     guess += nrhsentries;
891 :     }
892 :     if ( Rhstype[2] == 'X' ) {
893 :     for ( j=0;j<nrhsentries;j++ ) {
894 :     fprintf(out_file,rformat,exact[j]);
895 :     if ( acount++%Rhsperline == linemod ) fprintf(out_file,"\n");
896 :     }
897 :     if ( acount%Rhsperline != linemod ) {
898 :     fprintf(out_file,"\n");
899 :     linemod = (acount-1)%Rhsperline;
900 :     }
901 :     exact += nrhsentries;
902 :     }
903 :     }
904 :     }
905 :    
906 :     }
907 :    
908 :     if ( fclose(out_file) != 0){
909 :     fprintf(stderr,"Error closing file in writeHB_mat_double().\n");
910 :     return 0;
911 :     } else return 1;
912 :    
913 :     }
914 :    
915 :     int readHB_mat_char(const char* filename, int colptr[], int rowind[],
916 :     char val[], char* Valfmt)
917 :     {
918 :     /****************************************************************************/
919 :     /* This function opens and reads the specified file, interpreting its */
920 :     /* contents as a sparse matrix stored in the Harwell/Boeing standard */
921 :     /* format and creating compressed column storage scheme vectors to hold */
922 :     /* the index and nonzero value information. */
923 :     /* */
924 :     /* ---------- */
925 :     /* **CAVEAT** */
926 :     /* ---------- */
927 :     /* Parsing real formats from Fortran is tricky, and this file reader */
928 :     /* does not claim to be foolproof. It has been tested for cases when */
929 :     /* the real values are printed consistently and evenly spaced on each */
930 :     /* line, with Fixed (F), and Exponential (E or D) formats. */
931 :     /* */
932 :     /* ** If the input file does not adhere to the H/B format, the ** */
933 :     /* ** results will be unpredictable. ** */
934 :     /* */
935 :     /****************************************************************************/
936 :     FILE *in_file;
937 :     int i,j,ind,col,offset,count,last;
938 :     int Nrow,Ncol,Nnzero,Nentries,Nrhs;
939 :     int Ptrcrd, Indcrd, Valcrd, Rhscrd;
940 :     int Ptrperline, Ptrwidth, Indperline, Indwidth;
941 :     int Valperline, Valwidth, Valprec;
942 :     int Valflag; /* Indicates 'E','D', or 'F' float format */
943 :     char* ThisElement;
944 :     char line[BUFSIZ];
945 :     char Title[73], Key[8], Type[4], Rhstype[4];
946 :     char Ptrfmt[17], Indfmt[17], Rhsfmt[21];
947 :    
948 :     if ( (in_file = fopen( filename, "r")) == NULL ) {
949 :     fprintf(stderr,"Error: Cannot open file: %s\n",filename);
950 :     return 0;
951 :     }
952 :    
953 :     readHB_header(in_file, Title, Key, Type, &Nrow, &Ncol, &Nnzero, &Nrhs,
954 :     Ptrfmt, Indfmt, Valfmt, Rhsfmt,
955 :     &Ptrcrd, &Indcrd, &Valcrd, &Rhscrd, Rhstype);
956 :    
957 :     /* Parse the array input formats from Line 3 of HB file */
958 :     ParseIfmt(Ptrfmt,&Ptrperline,&Ptrwidth);
959 :     ParseIfmt(Indfmt,&Indperline,&Indwidth);
960 :     if ( Type[0] != 'P' ) { /* Skip if pattern only */
961 :     ParseRfmt(Valfmt,&Valperline,&Valwidth,&Valprec,&Valflag);
962 :     if (Valflag == 'D') {
963 :     *strchr(Valfmt,'D') = 'E';
964 :     }
965 :     }
966 :    
967 :     /* Read column pointer array: */
968 :    
969 :     offset = 1-_SP_base; /* if base 0 storage is declared (via macro definition), */
970 :     /* then storage entries are offset by 1 */
971 :    
972 :     ThisElement = (char *) malloc(Ptrwidth+1);
973 :     if ( ThisElement == NULL ) IOHBTerminate("Insufficient memory for ThisElement.");
974 :     *(ThisElement+Ptrwidth) = (char) NULL;
975 :     count=0;
976 :     for (i=0;i<Ptrcrd;i++)
977 :     {
978 :     fgets(line, BUFSIZ, in_file);
979 :     if ( sscanf(line,"%*s") < 0 )
980 :     IOHBTerminate("iohb.c: Null (or blank) line in pointer data region of HB file.\n");
981 :     col = 0;
982 :     for (ind = 0;ind<Ptrperline;ind++)
983 :     {
984 :     if (count > Ncol) break;
985 :     strncpy(ThisElement,line+col,Ptrwidth);
986 :     /*ThisElement = substr(line,col,Ptrwidth);*/
987 :     colptr[count] = atoi(ThisElement)-offset;
988 :     count++; col += Ptrwidth;
989 :     }
990 :     }
991 :     free(ThisElement);
992 :    
993 :     /* Read row index array: */
994 :    
995 :     ThisElement = (char *) malloc(Indwidth+1);
996 :     if ( ThisElement == NULL ) IOHBTerminate("Insufficient memory for ThisElement.");
997 :     *(ThisElement+Indwidth) = (char) NULL;
998 :     count = 0;
999 :     for (i=0;i<Indcrd;i++)
1000 :     {
1001 :     fgets(line, BUFSIZ, in_file);
1002 :     if ( sscanf(line,"%*s") < 0 )
1003 :     IOHBTerminate("iohb.c: Null (or blank) line in index data region of HB file.\n");
1004 :     col = 0;
1005 :     for (ind = 0;ind<Indperline;ind++)
1006 :     {
1007 :     if (count == Nnzero) break;
1008 :     strncpy(ThisElement,line+col,Indwidth);
1009 :     /*ThisElement = substr(line,col,Indwidth);*/
1010 :     rowind[count] = atoi(ThisElement)-offset;
1011 :     count++; col += Indwidth;
1012 :     }
1013 :     }
1014 :     free(ThisElement);
1015 :    
1016 :     /* Read array of values: AS CHARACTERS*/
1017 :    
1018 :     if ( Type[0] != 'P' ) { /* Skip if pattern only */
1019 :    
1020 :     if ( Type[0] == 'C' ) Nentries = 2*Nnzero;
1021 :     else Nentries = Nnzero;
1022 :    
1023 :     ThisElement = (char *) malloc(Valwidth+1);
1024 :     if ( ThisElement == NULL ) IOHBTerminate("Insufficient memory for ThisElement.");
1025 :     *(ThisElement+Valwidth) = (char) NULL;
1026 :     count = 0;
1027 :     for (i=0;i<Valcrd;i++)
1028 :     {
1029 :     fgets(line, BUFSIZ, in_file);
1030 :     if ( sscanf(line,"%*s") < 0 )
1031 :     IOHBTerminate("iohb.c: Null (or blank) line in value data region of HB file.\n");
1032 :     if (Valflag == 'D') {
1033 :     while( strchr(line,'D') ) *strchr(line,'D') = 'E';
1034 :     }
1035 :     col = 0;
1036 :     for (ind = 0;ind<Valperline;ind++)
1037 :     {
1038 :     if (count == Nentries) break;
1039 :     ThisElement = &val[count*Valwidth];
1040 :     strncpy(ThisElement,line+col,Valwidth);
1041 :     /*strncpy(ThisElement,substr(line,col,Valwidth),Valwidth);*/
1042 :     if ( Valflag != 'F' && strchr(ThisElement,'E') == NULL ) {
1043 :     /* insert a char prefix for exp */
1044 :     last = strlen(ThisElement);
1045 :     for (j=last+1;j>=0;j--) {
1046 :     ThisElement[j] = ThisElement[j-1];
1047 :     if ( ThisElement[j] == '+' || ThisElement[j] == '-' ) {
1048 :     ThisElement[j-1] = Valflag;
1049 :     break;
1050 :     }
1051 :     }
1052 :     }
1053 :     count++; col += Valwidth;
1054 :     }
1055 :     }
1056 :     }
1057 :    
1058 :     return 1;
1059 :     }
1060 :    
1061 :     int readHB_newmat_char(const char* filename, int* M, int* N, int* nonzeros, int** colptr,
1062 :     int** rowind, char** val, char** Valfmt)
1063 :     {
1064 :     FILE *in_file;
1065 :     int Nrhs;
1066 :     int Ptrcrd, Indcrd, Valcrd, Rhscrd;
1067 :     int Valperline, Valwidth, Valprec;
1068 :     int Valflag; /* Indicates 'E','D', or 'F' float format */
1069 :     char Title[73], Key[9], Type[4], Rhstype[4];
1070 :     char Ptrfmt[17], Indfmt[17], Rhsfmt[21];
1071 :    
1072 :     if ((in_file = fopen( filename, "r")) == NULL) {
1073 :     fprintf(stderr,"Error: Cannot open file: %s\n",filename);
1074 :     return 0;
1075 :     }
1076 :    
1077 :     *Valfmt = (char *)malloc(21*sizeof(char));
1078 :     if ( *Valfmt == NULL ) IOHBTerminate("Insufficient memory for Valfmt.");
1079 :     readHB_header(in_file, Title, Key, Type, M, N, nonzeros, &Nrhs,
1080 :     Ptrfmt, Indfmt, (*Valfmt), Rhsfmt,
1081 :     &Ptrcrd, &Indcrd, &Valcrd, &Rhscrd, Rhstype);
1082 :     fclose(in_file);
1083 :     ParseRfmt(*Valfmt,&Valperline,&Valwidth,&Valprec,&Valflag);
1084 :    
1085 :     *colptr = (int *)malloc((*N+1)*sizeof(int));
1086 :     if ( *colptr == NULL ) IOHBTerminate("Insufficient memory for colptr.\n");
1087 :     *rowind = (int *)malloc(*nonzeros*sizeof(int));
1088 :     if ( *rowind == NULL ) IOHBTerminate("Insufficient memory for rowind.\n");
1089 :     if ( Type[0] == 'C' ) {
1090 :     /*
1091 :     fprintf(stderr, "Warning: Reading complex data from HB file %s.\n",filename);
1092 :     fprintf(stderr, " Real and imaginary parts will be interlaced in val[].\n");
1093 :     */
1094 :     /* Malloc enough space for real AND imaginary parts of val[] */
1095 :     *val = (char *)malloc(*nonzeros*Valwidth*sizeof(char)*2);
1096 :     if ( *val == NULL ) IOHBTerminate("Insufficient memory for val.\n");
1097 :     } else {
1098 :     if ( Type[0] != 'P' ) {
1099 :     /* Malloc enough space for real array val[] */
1100 :     *val = (char *)malloc(*nonzeros*Valwidth*sizeof(char));
1101 :     if ( *val == NULL ) IOHBTerminate("Insufficient memory for val.\n");
1102 :     }
1103 :     } /* No val[] space needed if pattern only */
1104 :     return readHB_mat_char(filename, *colptr, *rowind, *val, *Valfmt);
1105 :    
1106 :     }
1107 :    
1108 :     int readHB_aux_char(const char* filename, const char AuxType, char b[])
1109 :     {
1110 :     /****************************************************************************/
1111 :     /* This function opens and reads the specified file, placing auxilary */
1112 :     /* vector(s) of the given type (if available) in b : */
1113 :     /* Return value is the number of vectors successfully read. */
1114 :     /* */
1115 :     /* AuxType = 'F' full right-hand-side vector(s) */
1116 :     /* AuxType = 'G' initial Guess vector(s) */
1117 :     /* AuxType = 'X' eXact solution vector(s) */
1118 :     /* */
1119 :     /* ---------- */
1120 :     /* **CAVEAT** */
1121 :     /* ---------- */
1122 :     /* Parsing real formats from Fortran is tricky, and this file reader */
1123 :     /* does not claim to be foolproof. It has been tested for cases when */
1124 :     /* the real values are printed consistently and evenly spaced on each */
1125 :     /* line, with Fixed (F), and Exponential (E or D) formats. */
1126 :     /* */
1127 :     /* ** If the input file does not adhere to the H/B format, the ** */
1128 :     /* ** results will be unpredictable. ** */
1129 :     /* */
1130 :     /****************************************************************************/
1131 :     FILE *in_file;
1132 :     int i,j,n,maxcol,start,stride,col,last,linel,nvecs,rhsi;
1133 :     int Nrow, Ncol, Nnzero, Nentries,Nrhs;
1134 :     int Ptrcrd, Indcrd, Valcrd, Rhscrd;
1135 :     int Rhsperline, Rhswidth, Rhsprec;
1136 :     int Rhsflag;
1137 :     char Title[73], Key[9], Type[4], Rhstype[4];
1138 :     char Ptrfmt[17], Indfmt[17], Valfmt[21], Rhsfmt[21];
1139 :     char line[BUFSIZ];
1140 :     char *ThisElement;
1141 :    
1142 :     if ((in_file = fopen( filename, "r")) == NULL) {
1143 :     fprintf(stderr,"Error: Cannot open file: %s\n",filename);
1144 :     return 0;
1145 :     }
1146 :    
1147 :     readHB_header(in_file, Title, Key, Type, &Nrow, &Ncol, &Nnzero, &Nrhs,
1148 :     Ptrfmt, Indfmt, Valfmt, Rhsfmt,
1149 :     &Ptrcrd, &Indcrd, &Valcrd, &Rhscrd, Rhstype);
1150 :    
1151 :     if (Nrhs <= 0)
1152 :     {
1153 :     fprintf(stderr, "Warn: Attempt to read auxillary vector(s) when none are present.\n");
1154 :     return 0;
1155 :     }
1156 :     if (Rhstype[0] != 'F' )
1157 :     {
1158 :     fprintf(stderr,"Warn: Attempt to read auxillary vector(s) which are not stored in Full form.\n");
1159 :     fprintf(stderr," Rhs must be specified as full. \n");
1160 :     return 0;
1161 :     }
1162 :    
1163 :     /* If reading complex data, allow for interleaved real and imaginary values. */
1164 :     if ( Type[0] == 'C' ) {
1165 :     Nentries = 2*Nrow;
1166 :     } else {
1167 :     Nentries = Nrow;
1168 :     }
1169 :    
1170 :     nvecs = 1;
1171 :    
1172 :     if ( Rhstype[1] == 'G' ) nvecs++;
1173 :     if ( Rhstype[2] == 'X' ) nvecs++;
1174 :    
1175 :     if ( AuxType == 'G' && Rhstype[1] != 'G' ) {
1176 :     fprintf(stderr, "Warn: Attempt to read auxillary Guess vector(s) when none are present.\n");
1177 :     return 0;
1178 :     }
1179 :     if ( AuxType == 'X' && Rhstype[2] != 'X' ) {
1180 :     fprintf(stderr, "Warn: Attempt to read auxillary eXact solution vector(s) when none are present.\n");
1181 :     return 0;
1182 :     }
1183 :    
1184 :     ParseRfmt(Rhsfmt, &Rhsperline, &Rhswidth, &Rhsprec,&Rhsflag);
1185 :     maxcol = Rhsperline*Rhswidth;
1186 :    
1187 :     /* Lines to skip before starting to read RHS values... */
1188 :     n = Ptrcrd + Indcrd + Valcrd;
1189 :    
1190 :     for (i = 0; i < n; i++)
1191 :     fgets(line, BUFSIZ, in_file);
1192 :    
1193 :     /* start - number of initial aux vector entries to skip */
1194 :     /* to reach first vector requested */
1195 :     /* stride - number of aux vector entries to skip between */
1196 :     /* requested vectors */
1197 :     if ( AuxType == 'F' ) start = 0;
1198 :     else if ( AuxType == 'G' ) start = Nentries;
1199 :     else start = (nvecs-1)*Nentries;
1200 :     stride = (nvecs-1)*Nentries;
1201 :    
1202 :     fgets(line, BUFSIZ, in_file);
1203 :     linel= strchr(line,'\n')-line;
1204 :     if ( sscanf(line,"%*s") < 0 )
1205 :     IOHBTerminate("iohb.c: Null (or blank) line in auxillary vector data region of HB file.\n");
1206 :     col = 0;
1207 :     /* Skip to initial offset */
1208 :    
1209 :     for (i=0;i<start;i++) {
1210 :     col += Rhswidth;
1211 :     if ( col >= ( maxcol<linel?maxcol:linel ) ) {
1212 :     fgets(line, BUFSIZ, in_file);
1213 :     linel= strchr(line,'\n')-line;
1214 :     if ( sscanf(line,"%*s") < 0 )
1215 :     IOHBTerminate("iohb.c: Null (or blank) line in auxillary vector data region of HB file.\n");
1216 :     col = 0;
1217 :     }
1218 :     }
1219 :    
1220 :     if (Rhsflag == 'D') {
1221 :     while( strchr(line,'D') ) *strchr(line,'D') = 'E';
1222 :     }
1223 :     /* Read a vector of desired type, then skip to next */
1224 :     /* repeating to fill Nrhs vectors */
1225 :    
1226 :     for (rhsi=0;rhsi<Nrhs;rhsi++) {
1227 :    
1228 :     for (i=0;i<Nentries;i++) {
1229 :     if ( col >= ( maxcol<linel?maxcol:linel ) ) {
1230 :     fgets(line, BUFSIZ, in_file);
1231 :     linel= strchr(line,'\n')-line;
1232 :     if ( sscanf(line,"%*s") < 0 )
1233 :     IOHBTerminate("iohb.c: Null (or blank) line in auxillary vector data region of HB file.\n");
1234 :     if (Rhsflag == 'D') {
1235 :     while( strchr(line,'D') ) *strchr(line,'D') = 'E';
1236 :     }
1237 :     col = 0;
1238 :     }
1239 :     ThisElement = &b[i*Rhswidth];
1240 :     strncpy(ThisElement,line+col,Rhswidth);
1241 :     if ( Rhsflag != 'F' && strchr(ThisElement,'E') == NULL ) {
1242 :     /* insert a char prefix for exp */
1243 :     last = strlen(ThisElement);
1244 :     for (j=last+1;j>=0;j--) {
1245 :     ThisElement[j] = ThisElement[j-1];
1246 :     if ( ThisElement[j] == '+' || ThisElement[j] == '-' ) {
1247 :     ThisElement[j-1] = Rhsflag;
1248 :     break;
1249 :     }
1250 :     }
1251 :     }
1252 :     col += Rhswidth;
1253 :     }
1254 :     b+=Nentries*Rhswidth;
1255 :    
1256 :     /* Skip any interleaved Guess/eXact vectors */
1257 :    
1258 :     for (i=0;i<stride;i++) {
1259 :     col += Rhswidth;
1260 :     if ( col >= ( maxcol<linel?maxcol:linel ) ) {
1261 :     fgets(line, BUFSIZ, in_file);
1262 :     linel= strchr(line,'\n')-line;
1263 :     if ( sscanf(line,"%*s") < 0 )
1264 :     IOHBTerminate("iohb.c: Null (or blank) line in auxillary vector data region of HB file.\n");
1265 :     col = 0;
1266 :     }
1267 :     }
1268 :    
1269 :     }
1270 :    
1271 :    
1272 :     fclose(in_file);
1273 :     return Nrhs;
1274 :     }
1275 :    
1276 :     int readHB_newaux_char(const char* filename, const char AuxType, char** b, char** Rhsfmt)
1277 :     {
1278 :     FILE *in_file;
1279 :     int Ptrcrd, Indcrd, Valcrd, Rhscrd;
1280 :     int Nrow,Ncol,Nnzero,Nrhs;
1281 :     int Rhsperline, Rhswidth, Rhsprec;
1282 :     int Rhsflag;
1283 :     char Title[73], Key[9], Type[4], Rhstype[4];
1284 :     char Ptrfmt[17], Indfmt[17], Valfmt[21];
1285 :    
1286 :     if ((in_file = fopen( filename, "r")) == NULL) {
1287 :     fprintf(stderr,"Error: Cannot open file: %s\n",filename);
1288 :     return 0;
1289 :     }
1290 :    
1291 :     *Rhsfmt = (char *)malloc(21*sizeof(char));
1292 :     if ( *Rhsfmt == NULL ) IOHBTerminate("Insufficient memory for Rhsfmt.");
1293 :     readHB_header(in_file, Title, Key, Type, &Nrow, &Ncol, &Nnzero, &Nrhs,
1294 :     Ptrfmt, Indfmt, Valfmt, (*Rhsfmt),
1295 :     &Ptrcrd, &Indcrd, &Valcrd, &Rhscrd, Rhstype);
1296 :     fclose(in_file);
1297 :     if ( Nrhs == 0 ) {
1298 :     fprintf(stderr,"Warn: Requested read of aux vector(s) when none are present.\n");
1299 :     return 0;
1300 :     } else {
1301 :     ParseRfmt(*Rhsfmt,&Rhsperline,&Rhswidth,&Rhsprec,&Rhsflag);
1302 :     if ( Type[0] == 'C' ) {
1303 :     fprintf(stderr, "Warning: Reading complex aux vector(s) from HB file %s.",filename);
1304 :     fprintf(stderr, " Real and imaginary parts will be interlaced in b[].");
1305 :     *b = (char *)malloc(Nrow*Nrhs*Rhswidth*sizeof(char)*2);
1306 :     if ( *b == NULL ) IOHBTerminate("Insufficient memory for rhs.\n");
1307 :     return readHB_aux_char(filename, AuxType, *b);
1308 :     } else {
1309 :     *b = (char *)malloc(Nrow*Nrhs*Rhswidth*sizeof(char));
1310 :     if ( *b == NULL ) IOHBTerminate("Insufficient memory for rhs.\n");
1311 :     return readHB_aux_char(filename, AuxType, *b);
1312 :     }
1313 :     }
1314 :     }
1315 :    
1316 :     int writeHB_mat_char(const char* filename, int M, int N,
1317 :     int nz, const int colptr[], const int rowind[],
1318 :     const char val[], int Nrhs, const char rhs[],
1319 :     const char guess[], const char exact[],
1320 :     const char* Title, const char* Key, const char* Type,
1321 :     char* Ptrfmt, char* Indfmt, char* Valfmt, char* Rhsfmt,
1322 :     const char* Rhstype)
1323 :     {
1324 :     /****************************************************************************/
1325 :     /* The writeHB function opens the named file and writes the specified */
1326 :     /* matrix and optional right-hand-side(s) to that file in Harwell-Boeing */
1327 :     /* format. */
1328 :     /* */
1329 :     /* For a description of the Harwell Boeing standard, see: */
1330 :     /* Duff, et al., ACM TOMS Vol.15, No.1, March 1989 */
1331 :     /* */
1332 :     /****************************************************************************/
1333 :     FILE *out_file;
1334 :     int i,j,acount,linemod,entry,offset;
1335 :     int totcrd, ptrcrd, indcrd, valcrd, rhscrd;
1336 :     int nvalentries, nrhsentries;
1337 :     int Ptrperline, Ptrwidth, Indperline, Indwidth;
1338 :     int Rhsperline, Rhswidth, Rhsprec;
1339 :     int Rhsflag;
1340 :     int Valperline, Valwidth, Valprec;
1341 :     int Valflag; /* Indicates 'E','D', or 'F' float format */
1342 :     char pformat[16],iformat[16],vformat[19],rformat[19];
1343 :    
1344 :     if ( Type[0] == 'C' ) {
1345 :     nvalentries = 2*nz;
1346 :     nrhsentries = 2*M;
1347 :     } else {
1348 :     nvalentries = nz;
1349 :     nrhsentries = M;
1350 :     }
1351 :    
1352 :     if ( filename != NULL ) {
1353 :     if ( (out_file = fopen( filename, "w")) == NULL ) {
1354 :     fprintf(stderr,"Error: Cannot open file: %s\n",filename);
1355 :     return 0;
1356 :     }
1357 :     } else out_file = stdout;
1358 :    
1359 :     if ( Ptrfmt == NULL ) Ptrfmt = "(8I10)";
1360 :     ParseIfmt(Ptrfmt,&Ptrperline,&Ptrwidth);
1361 :     sprintf(pformat,"%%%dd",Ptrwidth);
1362 :    
1363 :     if ( Indfmt == NULL ) Indfmt = Ptrfmt;
1364 :     ParseIfmt(Indfmt,&Indperline,&Indwidth);
1365 :     sprintf(iformat,"%%%dd",Indwidth);
1366 :    
1367 :     if ( Type[0] != 'P' ) { /* Skip if pattern only */
1368 :     if ( Valfmt == NULL ) Valfmt = "(4E20.13)";
1369 :     ParseRfmt(Valfmt,&Valperline,&Valwidth,&Valprec,&Valflag);
1370 :     sprintf(vformat,"%%%ds",Valwidth);
1371 :     }
1372 :    
1373 :     ptrcrd = (N+1)/Ptrperline;
1374 :     if ( (N+1)%Ptrperline != 0) ptrcrd++;
1375 :    
1376 :     indcrd = nz/Indperline;
1377 :     if ( nz%Indperline != 0) indcrd++;
1378 :    
1379 :     valcrd = nvalentries/Valperline;
1380 :     if ( nvalentries%Valperline != 0) valcrd++;
1381 :    
1382 :     if ( Nrhs > 0 ) {
1383 :     if ( Rhsfmt == NULL ) Rhsfmt = Valfmt;
1384 :     ParseRfmt(Rhsfmt,&Rhsperline,&Rhswidth,&Rhsprec, &Rhsflag);
1385 :     sprintf(rformat,"%%%ds",Rhswidth);
1386 :     rhscrd = nrhsentries/Rhsperline;
1387 :     if ( nrhsentries%Rhsperline != 0) rhscrd++;
1388 :     if ( Rhstype[1] == 'G' ) rhscrd+=rhscrd;
1389 :     if ( Rhstype[2] == 'X' ) rhscrd+=rhscrd;
1390 :     rhscrd*=Nrhs;
1391 :     } else rhscrd = 0;
1392 :    
1393 :     totcrd = 4+ptrcrd+indcrd+valcrd+rhscrd;
1394 :    
1395 :    
1396 :     /* Print header information: */
1397 :    
1398 :     fprintf(out_file,"%-72s%-8s\n%14d%14d%14d%14d%14d\n",Title, Key, totcrd,
1399 :     ptrcrd, indcrd, valcrd, rhscrd);
1400 :     fprintf(out_file,"%3s%11s%14d%14d%14d\n",Type," ", M, N, nz);
1401 :     fprintf(out_file,"%-16s%-16s%-20s", Ptrfmt, Indfmt, Valfmt);
1402 :     if ( Nrhs != 0 ) {
1403 :     /* Print Rhsfmt on fourth line and */
1404 :     /* optional fifth header line for auxillary vector information: */
1405 :     fprintf(out_file,"%-20s\n%-14s%d\n",Rhsfmt,Rhstype,Nrhs);
1406 :     } else fprintf(out_file,"\n");
1407 :    
1408 :     offset = 1-_SP_base; /* if base 0 storage is declared (via macro definition), */
1409 :     /* then storage entries are offset by 1 */
1410 :    
1411 :     /* Print column pointers: */
1412 :     for (i=0;i<N+1;i++)
1413 :     {
1414 :     entry = colptr[i]+offset;
1415 :     fprintf(out_file,pformat,entry);
1416 :     if ( (i+1)%Ptrperline == 0 ) fprintf(out_file,"\n");
1417 :     }
1418 :    
1419 :     if ( (N+1) % Ptrperline != 0 ) fprintf(out_file,"\n");
1420 :    
1421 :     /* Print row indices: */
1422 :     for (i=0;i<nz;i++)
1423 :     {
1424 :     entry = rowind[i]+offset;
1425 :     fprintf(out_file,iformat,entry);
1426 :     if ( (i+1)%Indperline == 0 ) fprintf(out_file,"\n");
1427 :     }
1428 :    
1429 :     if ( nz % Indperline != 0 ) fprintf(out_file,"\n");
1430 :    
1431 :     /* Print values: */
1432 :    
1433 :     if ( Type[0] != 'P' ) { /* Skip if pattern only */
1434 :     for (i=0;i<nvalentries;i++)
1435 :     {
1436 :     fprintf(out_file,vformat,val+i*Valwidth);
1437 :     if ( (i+1)%Valperline == 0 ) fprintf(out_file,"\n");
1438 :     }
1439 :    
1440 :     if ( nvalentries % Valperline != 0 ) fprintf(out_file,"\n");
1441 :    
1442 :     /* Print right hand sides: */
1443 :     acount = 1;
1444 :     linemod=0;
1445 :     if ( Nrhs > 0 ) {
1446 :     for (j=0;j<Nrhs;j++) {
1447 :     for (i=0;i<nrhsentries;i++)
1448 :     {
1449 :     fprintf(out_file,rformat,rhs+i*Rhswidth);
1450 :     if ( acount++%Rhsperline == linemod ) fprintf(out_file,"\n");
1451 :     }
1452 :     if ( acount%Rhsperline != linemod ) {
1453 :     fprintf(out_file,"\n");
1454 :     linemod = (acount-1)%Rhsperline;
1455 :     }
1456 :     if ( Rhstype[1] == 'G' ) {
1457 :     for (i=0;i<nrhsentries;i++)
1458 :     {
1459 :     fprintf(out_file,rformat,guess+i*Rhswidth);
1460 :     if ( acount++%Rhsperline == linemod ) fprintf(out_file,"\n");
1461 :     }
1462 :     if ( acount%Rhsperline != linemod ) {
1463 :     fprintf(out_file,"\n");
1464 :     linemod = (acount-1)%Rhsperline;
1465 :     }
1466 :     }
1467 :     if ( Rhstype[2] == 'X' ) {
1468 :     for (i=0;i<nrhsentries;i++)
1469 :     {
1470 :     fprintf(out_file,rformat,exact+i*Rhswidth);
1471 :     if ( acount++%Rhsperline == linemod ) fprintf(out_file,"\n");
1472 :     }
1473 :     if ( acount%Rhsperline != linemod ) {
1474 :     fprintf(out_file,"\n");
1475 :     linemod = (acount-1)%Rhsperline;
1476 :     }
1477 :     }
1478 :     }
1479 :     }
1480 :    
1481 :     }
1482 :    
1483 :     if ( fclose(out_file) != 0){
1484 :     fprintf(stderr,"Error closing file in writeHB_mat_char().\n");
1485 :     return 0;
1486 :     } else return 1;
1487 :    
1488 :     }
1489 :    
1490 :     int ParseIfmt(char* fmt, int* perline, int* width)
1491 :     {
1492 :     /*************************************************/
1493 :     /* Parse an *integer* format field to determine */
1494 :     /* width and number of elements per line. */
1495 :     /*************************************************/
1496 :     char *tmp;
1497 :     if (fmt == NULL ) {
1498 :     *perline = 0; *width = 0; return 0;
1499 :     }
1500 :     upcase(fmt);
1501 :     tmp = strchr(fmt,'(');
1502 :     tmp = substr(fmt,tmp - fmt + 1, strchr(fmt,'I') - tmp - 1);
1503 :     *perline = atoi(tmp);
1504 :     tmp = strchr(fmt,'I');
1505 :     tmp = substr(fmt,tmp - fmt + 1, strchr(fmt,')') - tmp - 1);
1506 :     return *width = atoi(tmp);
1507 :     }
1508 :    
1509 :     int ParseRfmt(char* fmt, int* perline, int* width, int* prec, int* flag)
1510 :     {
1511 :     /*************************************************/
1512 :     /* Parse a *real* format field to determine */
1513 :     /* width and number of elements per line. */
1514 :     /* Also sets flag indicating 'E' 'F' 'P' or 'D' */
1515 :     /* format. */
1516 :     /*************************************************/
1517 :     char* tmp;
1518 :     char* tmp2;
1519 :     char* tmp3;
1520 :     int len;
1521 :    
1522 :     if (fmt == NULL ) {
1523 :     *perline = 0;
1524 :     *width = 0;
1525 :     flag = NULL;
1526 :     return 0;
1527 :     }
1528 :    
1529 :     upcase(fmt);
1530 :     if (strchr(fmt,'(') != NULL) fmt = strchr(fmt,'(');
1531 :     if (strchr(fmt,')') != NULL) {
1532 :     tmp2 = strchr(fmt,')');
1533 :     while ( strchr(tmp2+1,')') != NULL ) {
1534 :     tmp2 = strchr(tmp2+1,')');
1535 :     }
1536 :     *(tmp2+1) = (int) NULL;
1537 :     }
1538 :     if (strchr(fmt,'P') != NULL) /* Remove any scaling factor, which */
1539 :     { /* affects output only, not input */
1540 :     if (strchr(fmt,'(') != NULL) {
1541 :     tmp = strchr(fmt,'P');
1542 :     if ( *(++tmp) == ',' ) tmp++;
1543 :     tmp3 = strchr(fmt,'(')+1;
1544 :     len = tmp-tmp3;
1545 :     tmp2 = tmp3;
1546 :     while ( *(tmp2+len) != (int) NULL ) {
1547 :     *tmp2=*(tmp2+len);
1548 :     tmp2++;
1549 :     }
1550 :     *(strchr(fmt,')')+1) = (int) NULL;
1551 :     }
1552 :     }
1553 :     if (strchr(fmt,'E') != NULL) {
1554 :     *flag = 'E';
1555 :     } else if (strchr(fmt,'D') != NULL) {
1556 :     *flag = 'D';
1557 :     } else if (strchr(fmt,'F') != NULL) {
1558 :     *flag = 'F';
1559 :     } else {
1560 :     fprintf(stderr,"Real format %s in H/B file not supported.\n",fmt);
1561 :     return 0;
1562 :     }
1563 :     tmp = strchr(fmt,'(');
1564 :     tmp = substr(fmt,tmp - fmt + 1, strchr(fmt,*flag) - tmp - 1);
1565 :     *perline = atoi(tmp);
1566 :     tmp = strchr(fmt,*flag);
1567 :     if ( strchr(fmt,'.') ) {
1568 :     *prec = atoi( substr( fmt, strchr(fmt,'.') - fmt + 1, strchr(fmt,')') - strchr(fmt,'.')-1) );
1569 :     tmp = substr(fmt,tmp - fmt + 1, strchr(fmt,'.') - tmp - 1);
1570 :     } else {
1571 :     tmp = substr(fmt,tmp - fmt + 1, strchr(fmt,')') - tmp - 1);
1572 :     }
1573 :     return *width = atoi(tmp);
1574 :     }
1575 :    
1576 :     char* substr(const char* S, const int pos, const int len)
1577 :     {
1578 :     int i;
1579 :     char *SubS;
1580 :     if ( pos+len <= strlen(S)) {
1581 :     SubS = (char *)malloc(len+1);
1582 :     if ( SubS == NULL ) IOHBTerminate("Insufficient memory for SubS.");
1583 :     for (i=0;i<len;i++) SubS[i] = S[pos+i];
1584 :     SubS[len] = (char) NULL;
1585 :     } else {
1586 :     SubS = NULL;
1587 :     }
1588 :     return SubS;
1589 :     }
1590 :    
1591 :     #include<ctype.h>
1592 :     void upcase(char* S)
1593 :     {
1594 :     /* Convert S to uppercase */
1595 :     int i,len;
1596 :     len = strlen(S);
1597 :     for (i=0;i< len;i++)
1598 :     S[i] = toupper(S[i]);
1599 :     }
1600 :    
1601 :     void IOHBTerminate(char* message)
1602 :     {
1603 :     fprintf(stderr,message);
1604 :     exit(1);
1605 :     }
1606 :    

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