Tpetra parallel linear algebra  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Tpetra_Util_iohb.cpp
1 /*
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(std::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 #ifndef _SP_base
213 #define _SP_base 1
214 #endif
215 /*---------------------------------------------------------------------*/
216 
217 #include "Tpetra_Util_iohb.h"
218 
219 #include <cstring>
220 #include <cmath>
221 #include <cstdlib>
222 using std::free;
223 using std::malloc;
224 using std::size_t;
225 
226 char* substr(const char* S, const int pos, const int len);
227 void upcase(char* S);
228 void IOHBTerminate(const char* message);
229 
230 int readHB_info(const char* filename, int* M, int* N, int* nz, char** Type,
231  int* Nrhs) {
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  std::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_type\n");
258 
259  if ((in_file = std::fopen(filename, "r")) == NULL) {
260  std::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  std::fclose(in_file);
268  *Type = mat_type;
269  *(*Type + 3) = '\0';
270  *M = Nrow;
271  *N = Ncol;
272  *nz = Nnzero;
273  if (Rhscrd == 0) {
274  *Nrhs = 0;
275  }
276 
277  /* In verbose mode, print some of the header information: */
278  /*
279  if (verbose == 1)
280  {
281  printf("Reading from Harwell-Boeing file %s (verbose on)...\n",filename);
282  printf(" Title: %s\n",Title);
283  printf(" Key: %s\n",Key);
284  printf(" The stored matrix is %i by %i with %i nonzeros.\n",
285  *M, *N, *nz );
286  printf(" %i right-hand--side(s) stored.\n",*Nrhs);
287  }
288  */
289 
290  return 1;
291 }
292 
293 int readHB_header(std::FILE* in_file, char* Title, char* Key, char* Type,
294  int* Nrow, int* Ncol, int* Nnzero, int* Nrhs,
295  char* Ptrfmt, char* Indfmt, char* Valfmt, char* Rhsfmt,
296  int* Ptrcrd, int* Indcrd, int* Valcrd, int* Rhscrd,
297  char* Rhstype) {
298  /*************************************************************************/
299  /* Read header information from the named H/B file... */
300  /*************************************************************************/
301  int Totcrd, Neltvl, Nrhsix;
302  char line[BUFSIZ];
303 
304  /* First line: */
305  if (std::fgets(line, BUFSIZ, in_file) == NULL) {
306  std::fprintf(stderr, "Error: Failed to read from file.\n");
307  return 0;
308  }
309  if (std::sscanf(line, "%*s") < 0)
310  IOHBTerminate("Trilinos_Util_iohb.cpp: Null (or blank) first line of HB file.\n");
311  (void)std::sscanf(line, "%72c%8[^\n]", Title, Key);
312  *(Key + 8) = '\0';
313  *(Title + 72) = '\0';
314 
315  /* Second line: */
316  if (std::fgets(line, BUFSIZ, in_file) == NULL) {
317  std::fprintf(stderr, "Error: Failed to read from file.\n");
318  return 0;
319  }
320  if (std::sscanf(line, "%*s") < 0)
321  IOHBTerminate("Trilinos_Util_iohb.cpp: Null (or blank) second line of HB file.\n");
322  if (std::sscanf(line, "%i", &Totcrd) != 1) Totcrd = 0;
323  if (std::sscanf(line, "%*i%i", Ptrcrd) != 1) *Ptrcrd = 0;
324  if (std::sscanf(line, "%*i%*i%i", Indcrd) != 1) *Indcrd = 0;
325  if (std::sscanf(line, "%*i%*i%*i%i", Valcrd) != 1) *Valcrd = 0;
326  if (std::sscanf(line, "%*i%*i%*i%*i%i", Rhscrd) != 1) *Rhscrd = 0;
327 
328  /* Third line: */
329  if (std::fgets(line, BUFSIZ, in_file) == NULL) {
330  std::fprintf(stderr, "Error: Failed to read from file.\n");
331  return 0;
332  }
333  if (std::sscanf(line, "%*s") < 0)
334  IOHBTerminate("Trilinos_Util_iohb.cpp: Null (or blank) third line of HB file.\n");
335  if (std::sscanf(line, "%3c", Type) != 1)
336  IOHBTerminate("Trilinos_Util_iohb.cpp: Invalid Type info, line 3 of Harwell-Boeing file.\n");
337  upcase(Type);
338  if (std::sscanf(line, "%*3c%i", Nrow) != 1) *Nrow = 0;
339  if (std::sscanf(line, "%*3c%*i%i", Ncol) != 1) *Ncol = 0;
340  if (std::sscanf(line, "%*3c%*i%*i%i", Nnzero) != 1) *Nnzero = 0;
341  if (std::sscanf(line, "%*3c%*i%*i%*i%i", &Neltvl) != 1) Neltvl = 0;
342 
343  /* Fourth line: */
344  if (std::fgets(line, BUFSIZ, in_file) == NULL) {
345  std::fprintf(stderr, "Error: Failed to read from file.\n");
346  return 0;
347  }
348  if (std::sscanf(line, "%*s") < 0)
349  IOHBTerminate("Trilinos_Util_iohb.cpp: Null (or blank) fourth line of HB file.\n");
350  if (std::sscanf(line, "%16c", Ptrfmt) != 1)
351  IOHBTerminate("Trilinos_Util_iohb.cpp: Invalid format info, line 4 of Harwell-Boeing file.\n");
352  if (std::sscanf(line, "%*16c%16c", Indfmt) != 1)
353  IOHBTerminate("Trilinos_Util_iohb.cpp: Invalid format info, line 4 of Harwell-Boeing file.\n");
354  if (std::sscanf(line, "%*16c%*16c%20c", Valfmt) != 1)
355  IOHBTerminate("Trilinos_Util_iohb.cpp: Invalid format info, line 4 of Harwell-Boeing file.\n");
356  std::sscanf(line, "%*16c%*16c%*20c%20c", Rhsfmt);
357  *(Ptrfmt + 16) = '\0';
358  *(Indfmt + 16) = '\0';
359  *(Valfmt + 20) = '\0';
360  *(Rhsfmt + 20) = '\0';
361 
362  /* (Optional) Fifth line: */
363  if (*Rhscrd != 0) {
364  if (std::fgets(line, BUFSIZ, in_file) == NULL) {
365  std::fprintf(stderr, "Error: Failed to read from file.\n");
366  return 0;
367  }
368  if (std::sscanf(line, "%*s") < 0)
369  IOHBTerminate("Trilinos_Util_iohb.cpp: Null (or blank) fifth line of HB file.\n");
370  if (std::sscanf(line, "%3c", Rhstype) != 1)
371  IOHBTerminate("Trilinos_Util_iohb.cpp: Invalid RHS type information, line 5 of Harwell-Boeing file.\n");
372  if (std::sscanf(line, "%*3c%i", Nrhs) != 1) *Nrhs = 0;
373  if (std::sscanf(line, "%*3c%*i%i", &Nrhsix) != 1) Nrhsix = 0;
374  }
375  return 1;
376 }
377 
378 int readHB_mat_double(const char* filename, int colptr[], int rowind[],
379  double val[]) {
380  /****************************************************************************/
381  /* This function opens and reads the specified file, interpreting its */
382  /* contents as a sparse matrix stored in the Harwell/Boeing standard */
383  /* format and creating compressed column storage scheme vectors to hold */
384  /* the index and nonzero value information. */
385  /* */
386  /* ---------- */
387  /* **CAVEAT** */
388  /* ---------- */
389  /* Parsing real formats from Fortran is tricky, and this file reader */
390  /* does not claim to be foolproof. It has been tested for cases when */
391  /* the real values are printed consistently and evenly spaced on each */
392  /* line, with Fixed (F), and Exponential (E or D) formats. */
393  /* */
394  /* ** If the input file does not adhere to the H/B format, the ** */
395  /* ** results will be unpredictable. ** */
396  /* */
397  /****************************************************************************/
398  std::FILE* in_file;
399  int i, j, ind, col, offset, count, last, Nrhs;
400  int Ptrcrd, Indcrd, Valcrd, Rhscrd;
401  int Nrow, Ncol, Nnzero, Nentries;
402  int Ptrperline, Ptrwidth, Indperline, Indwidth;
403  int Valperline, Valwidth, Valprec;
404  int Valflag; /* Indicates 'E','D', or 'F' float format */
405  char* ThisElement;
406  char Title[73], Key[9], Type[4] = "XXX", Rhstype[4];
407  char Ptrfmt[17], Indfmt[17], Valfmt[21], Rhsfmt[21];
408  char line[BUFSIZ];
409 
410  if ((in_file = std::fopen(filename, "r")) == NULL) {
411  std::fprintf(stderr, "Error: Cannot open file: %s\n", filename);
412  return 0;
413  }
414 
415  readHB_header(in_file, Title, Key, Type, &Nrow, &Ncol, &Nnzero, &Nrhs,
416  Ptrfmt, Indfmt, Valfmt, Rhsfmt,
417  &Ptrcrd, &Indcrd, &Valcrd, &Rhscrd, Rhstype);
418 
419  /* Parse the array input formats from Line 3 of HB file */
420  ParseIfmt(Ptrfmt, &Ptrperline, &Ptrwidth);
421  ParseIfmt(Indfmt, &Indperline, &Indwidth);
422  if (Type[0] != 'P') { /* Skip if pattern only */
423  ParseRfmt(Valfmt, &Valperline, &Valwidth, &Valprec, &Valflag);
424  }
425 
426  /* Read column pointer array: */
427 
428  offset = 1 - _SP_base; /* if base 0 storage is declared (via macro definition), */
429  /* then storage entries are offset by 1 */
430 
431  ThisElement = (char*)malloc(Ptrwidth + 1);
432  if (ThisElement == NULL) IOHBTerminate("Insufficient memory for ThisElement.");
433  *(ThisElement + Ptrwidth) = '\0';
434  count = 0;
435  for (i = 0; i < Ptrcrd; i++) {
436  if (std::fgets(line, BUFSIZ, in_file) == NULL) {
437  std::fprintf(stderr, "Error: Failed to read from file.\n");
438  return 0;
439  }
440  if (std::sscanf(line, "%*s") < 0)
441  IOHBTerminate("Trilinos_Util_iohb.cpp: Null (or blank) line in pointer data region of HB file.\n");
442  col = 0;
443  for (ind = 0; ind < Ptrperline; ind++) {
444  if (count > Ncol) break;
445  std::strncpy(ThisElement, line + col, Ptrwidth);
446  /* ThisElement = substr(line,col,Ptrwidth); */
447  colptr[count] = std::atoi(ThisElement) - offset;
448  count++;
449  col += Ptrwidth;
450  }
451  }
452  free(ThisElement);
453 
454  /* Read row index array: */
455 
456  ThisElement = (char*)malloc(Indwidth + 1);
457  if (ThisElement == NULL) IOHBTerminate("Insufficient memory for ThisElement.");
458  *(ThisElement + Indwidth) = '\0';
459  count = 0;
460  for (i = 0; i < Indcrd; i++) {
461  if (std::fgets(line, BUFSIZ, in_file) == NULL) {
462  std::fprintf(stderr, "Error: Failed to read from file.\n");
463  return 0;
464  }
465  if (std::sscanf(line, "%*s") < 0)
466  IOHBTerminate("Trilinos_Util_iohb.cpp: Null (or blank) line in index data region of HB file.\n");
467  col = 0;
468  for (ind = 0; ind < Indperline; ind++) {
469  if (count == Nnzero) break;
470  std::strncpy(ThisElement, line + col, Indwidth);
471  /* ThisElement = substr(line,col,Indwidth); */
472  rowind[count] = std::atoi(ThisElement) - offset;
473  count++;
474  col += Indwidth;
475  }
476  }
477  free(ThisElement);
478 
479  /* Read array of values: */
480 
481  if (Type[0] != 'P') { /* Skip if pattern only */
482 
483  if (Type[0] == 'C')
484  Nentries = 2 * Nnzero;
485  else
486  Nentries = Nnzero;
487 
488  ThisElement = (char*)malloc(Valwidth + 2);
489  if (ThisElement == NULL) IOHBTerminate("Insufficient memory for ThisElement.");
490  *(ThisElement + Valwidth) = '\0';
491  *(ThisElement + Valwidth + 1) = '\0';
492  count = 0;
493  for (i = 0; i < Valcrd; i++) {
494  if (std::fgets(line, BUFSIZ, in_file) == NULL) {
495  std::fprintf(stderr, "Error: Failed to read from file.\n");
496  return 0;
497  }
498  if (std::sscanf(line, "%*s") < 0)
499  IOHBTerminate("Trilinos_Util_iohb.cpp: Null (or blank) line in value data region of HB file.\n");
500  if (Valflag == 'D') {
501  while (std::strchr(line, 'D')) *std::strchr(line, 'D') = 'E';
502  /* *std::strchr(Valfmt,'D') = 'E'; */
503  }
504  col = 0;
505  for (ind = 0; ind < Valperline; ind++) {
506  if (count == Nentries) break;
507  std::strncpy(ThisElement, line + col, Valwidth);
508  /*ThisElement = substr(line,col,Valwidth);*/
509  if (Valflag != 'F' && std::strchr(ThisElement, 'E') == NULL) {
510  /* insert a char prefix for exp */
511  last = std::strlen(ThisElement);
512  for (j = last + 1; j >= 0; j--) {
513  ThisElement[j] = ThisElement[j - 1];
514  if (ThisElement[j] == '+' || ThisElement[j] == '-') {
515  ThisElement[j - 1] = Valflag;
516  break;
517  }
518  }
519  }
520  val[count] = std::atof(ThisElement);
521  count++;
522  col += Valwidth;
523  *(ThisElement + Valwidth) = '\0';
524  *(ThisElement + Valwidth + 1) = '\0';
525  }
526  }
527  free(ThisElement);
528  }
529 
530  std::fclose(in_file);
531  return 1;
532 }
533 
534 int readHB_newmat_double(const char* filename, int* M, int* N, int* nonzeros,
535  int** colptr, int** rowind, double** val) {
536  int Nrhs;
537  char* Type;
538 
539  if (readHB_info(filename, M, N, nonzeros, &Type, &Nrhs) == 0) {
540  return 0;
541  }
542 
543  *colptr = (int*)malloc((*N + 1) * sizeof(int));
544  if (*colptr == NULL) IOHBTerminate("Insufficient memory for colptr.\n");
545  *rowind = (int*)malloc(*nonzeros * sizeof(int));
546  if (*rowind == NULL) IOHBTerminate("Insufficient memory for rowind.\n");
547  if (Type[0] == 'C') {
548  /*
549  std::fprintf(stderr, "Warning: Reading complex data from HB file %s.\n",filename);
550  std::fprintf(stderr, " Real and imaginary parts will be interlaced in val[].\n");
551  */
552  /* Malloc enough space for real AND imaginary parts of val[] */
553  *val = (double*)malloc(*nonzeros * sizeof(double) * 2);
554  if (*val == NULL) IOHBTerminate("Insufficient memory for val.\n");
555  } else {
556  if (Type[0] != 'P') {
557  /* Malloc enough space for real array val[] */
558  *val = (double*)malloc(*nonzeros * sizeof(double));
559  if (*val == NULL) IOHBTerminate("Insufficient memory for val.\n");
560  }
561  } /* No val[] space needed if pattern only */
562  return readHB_mat_double(filename, *colptr, *rowind, *val);
563 }
564 
565 int readHB_aux_double(const char* filename, const char AuxType, double b[]) {
566  /****************************************************************************/
567  /* This function opens and reads the specified file, placing auxillary */
568  /* vector(s) of the given type (if available) in b. */
569  /* Return value is the number of vectors successfully read. */
570  /* */
571  /* AuxType = 'F' full right-hand-side vector(s) */
572  /* AuxType = 'G' initial Guess vector(s) */
573  /* AuxType = 'X' eXact solution vector(s) */
574  /* */
575  /* ---------- */
576  /* **CAVEAT** */
577  /* ---------- */
578  /* Parsing real formats from Fortran is tricky, and this file reader */
579  /* does not claim to be foolproof. It has been tested for cases when */
580  /* the real values are printed consistently and evenly spaced on each */
581  /* line, with Fixed (F), and Exponential (E or D) formats. */
582  /* */
583  /* ** If the input file does not adhere to the H/B format, the ** */
584  /* ** results will be unpredictable. ** */
585  /* */
586  /****************************************************************************/
587  std::FILE* in_file;
588  int i, j, n, maxcol, start, stride, col, last, linel;
589  int Ptrcrd, Indcrd, Valcrd, Rhscrd;
590  int Nrow, Ncol, Nnzero, Nentries;
591  int Nrhs, nvecs, rhsi;
592  int Rhsperline, Rhswidth, Rhsprec;
593  int Rhsflag;
594  char* ThisElement;
595  char Title[73], Key[9], Type[4] = "XXX", Rhstype[4];
596  char Ptrfmt[17], Indfmt[17], Valfmt[21], Rhsfmt[21];
597  char line[BUFSIZ];
598 
599  if ((in_file = std::fopen(filename, "r")) == NULL) {
600  std::fprintf(stderr, "Error: Cannot open file: %s\n", filename);
601  return 0;
602  }
603 
604  readHB_header(in_file, Title, Key, Type, &Nrow, &Ncol, &Nnzero, &Nrhs,
605  Ptrfmt, Indfmt, Valfmt, Rhsfmt,
606  &Ptrcrd, &Indcrd, &Valcrd, &Rhscrd, Rhstype);
607 
608  if (Nrhs <= 0) {
609  std::fprintf(stderr, "Warn: Attempt to read auxillary vector(s) when none are present.\n");
610  return 0;
611  }
612  if (Rhstype[0] != 'F') {
613  std::fprintf(stderr, "Warn: Attempt to read auxillary vector(s) which are not stored in Full form.\n");
614  std::fprintf(stderr, " Rhs must be specified as full. \n");
615  return 0;
616  }
617 
618  /* If reading complex data, allow for interleaved real and imaginary values. */
619  if (Type[0] == 'C') {
620  Nentries = 2 * Nrow;
621  } else {
622  Nentries = Nrow;
623  }
624 
625  nvecs = 1;
626 
627  if (Rhstype[1] == 'G') nvecs++;
628  if (Rhstype[2] == 'X') nvecs++;
629 
630  if (AuxType == 'G' && Rhstype[1] != 'G') {
631  std::fprintf(stderr, "Warn: Attempt to read auxillary Guess vector(s) when none are present.\n");
632  return 0;
633  }
634  if (AuxType == 'X' && Rhstype[2] != 'X') {
635  std::fprintf(stderr, "Warn: Attempt to read auxillary eXact solution vector(s) when none are present.\n");
636  return 0;
637  }
638 
639  ParseRfmt(Rhsfmt, &Rhsperline, &Rhswidth, &Rhsprec, &Rhsflag);
640  maxcol = Rhsperline * Rhswidth;
641 
642  /* Lines to skip before starting to read RHS values... */
643  n = Ptrcrd + Indcrd + Valcrd;
644 
645  for (i = 0; i < n; i++) {
646  if (std::fgets(line, BUFSIZ, in_file) == NULL) {
647  std::fprintf(stderr, "Error: Failed to read from file.\n");
648  return 0;
649  }
650  }
651 
652  /* start - number of initial aux vector entries to skip */
653  /* to reach first vector requested */
654  /* stride - number of aux vector entries to skip between */
655  /* requested vectors */
656  if (AuxType == 'F')
657  start = 0;
658  else if (AuxType == 'G')
659  start = Nentries;
660  else
661  start = (nvecs - 1) * Nentries;
662  stride = (nvecs - 1) * Nentries;
663 
664  if (std::fgets(line, BUFSIZ, in_file) == NULL) {
665  std::fprintf(stderr, "Error: Failed to read from file.\n");
666  return 0;
667  }
668  linel = std::strchr(line, '\n') - line;
669  col = 0;
670  /* Skip to initial offset */
671 
672  for (i = 0; i < start; i++) {
673  if (col >= (maxcol < linel ? maxcol : linel)) {
674  if (std::fgets(line, BUFSIZ, in_file) == NULL) {
675  std::fprintf(stderr, "Error: Failed to read from file.\n");
676  return 0;
677  }
678  linel = std::strchr(line, '\n') - line;
679  col = 0;
680  }
681  col += Rhswidth;
682  }
683  if (Rhsflag == 'D') {
684  while (std::strchr(line, 'D')) *std::strchr(line, 'D') = 'E';
685  }
686 
687  /* Read a vector of desired type, then skip to next */
688  /* repeating to fill Nrhs vectors */
689 
690  ThisElement = (char*)malloc(Rhswidth + 1);
691  if (ThisElement == NULL) IOHBTerminate("Insufficient memory for ThisElement.");
692  *(ThisElement + Rhswidth) = '\0';
693  for (rhsi = 0; rhsi < Nrhs; rhsi++) {
694  for (i = 0; i < Nentries; i++) {
695  if (col >= (maxcol < linel ? maxcol : linel)) {
696  if (std::fgets(line, BUFSIZ, in_file) == NULL) {
697  std::fprintf(stderr, "Error: Failed to read from file.\n");
698  return 0;
699  }
700  linel = std::strchr(line, '\n') - line;
701  if (Rhsflag == 'D') {
702  while (std::strchr(line, 'D')) *std::strchr(line, 'D') = 'E';
703  }
704  col = 0;
705  }
706  std::strncpy(ThisElement, line + col, Rhswidth);
707  /*ThisElement = substr(line, col, Rhswidth);*/
708  if (Rhsflag != 'F' && std::strchr(ThisElement, 'E') == NULL) {
709  /* insert a char prefix for exp */
710  last = std::strlen(ThisElement);
711  for (j = last + 1; j >= 0; j--) {
712  ThisElement[j] = ThisElement[j - 1];
713  if (ThisElement[j] == '+' || ThisElement[j] == '-') {
714  ThisElement[j - 1] = Rhsflag;
715  break;
716  }
717  }
718  }
719  b[i] = std::atof(ThisElement);
720  col += Rhswidth;
721  }
722 
723  /* Skip any interleaved Guess/eXact vectors */
724 
725  for (i = 0; i < stride; i++) {
726  if (col >= (maxcol < linel ? maxcol : linel)) {
727  if (std::fgets(line, BUFSIZ, in_file) == NULL) {
728  std::fprintf(stderr, "Error: Failed to read from file.\n");
729  return 0;
730  }
731  linel = std::strchr(line, '\n') - line;
732  col = 0;
733  }
734  col += Rhswidth;
735  }
736  }
737  free(ThisElement);
738 
739  std::fclose(in_file);
740  return Nrhs;
741 }
742 
743 int readHB_newaux_double(const char* filename, const char AuxType, double** b) {
744  int Nrhs = 0;
745  int M = 0;
746  int N = 0;
747  int nonzeros = 0;
748  char* Type = NULL;
749 
750  readHB_info(filename, &M, &N, &nonzeros, &Type, &Nrhs);
751  if (Nrhs <= 0) {
752  std::fprintf(stderr, "Warn: Requested read of aux vector(s) when none are present.\n");
753  return 0;
754  } else {
755  if (Type[0] == 'C') {
756  std::fprintf(stderr, "Warning: Reading complex aux vector(s) from HB file %s.", filename);
757  std::fprintf(stderr, " Real and imaginary parts will be interlaced in b[].");
758  *b = (double*)malloc(M * Nrhs * sizeof(double) * 2);
759  if (*b == NULL) IOHBTerminate("Insufficient memory for rhs.\n");
760  return readHB_aux_double(filename, AuxType, *b);
761  } else {
762  *b = (double*)malloc(M * Nrhs * sizeof(double));
763  if (*b == NULL) IOHBTerminate("Insufficient memory for rhs.\n");
764  return readHB_aux_double(filename, AuxType, *b);
765  }
766  }
767 }
768 
769 int writeHB_mat_double(const char* filename, int M, int N,
770  int nz, const int colptr[], const int rowind[],
771  const double val[], int Nrhs, const double rhs[],
772  const double guess[], const double exact[],
773  const char* Title, const char* Key, const char* Type,
774  char* Ptrfmt, char* Indfmt, char* Valfmt, char* Rhsfmt,
775  const char* Rhstype) {
776  /****************************************************************************/
777  /* The writeHB function opens the named file and writes the specified */
778  /* matrix and optional right-hand-side(s) to that file in Harwell-Boeing */
779  /* format. */
780  /* */
781  /* For a description of the Harwell Boeing standard, see: */
782  /* Duff, et al., ACM TOMS Vol.15, No.1, March 1989 */
783  /* */
784  /****************************************************************************/
785  std::FILE* out_file;
786  int i, j, entry, offset, acount, linemod;
787  int totcrd, ptrcrd, indcrd, valcrd, rhscrd;
788  int nvalentries, nrhsentries;
789  int Ptrperline, Ptrwidth, Indperline, Indwidth;
790  int Rhsperline, Rhswidth, Rhsprec;
791  int Rhsflag;
792  int Valperline, Valwidth, Valprec;
793  int Valflag; /* Indicates 'E','D', or 'F' float format */
794  char pformat[16], iformat[16], vformat[19], rformat[19];
795 
796  if (Type[0] == 'C') {
797  nvalentries = 2 * nz;
798  nrhsentries = 2 * M;
799  } else {
800  nvalentries = nz;
801  nrhsentries = M;
802  }
803 
804  if (filename != NULL) {
805  if ((out_file = std::fopen(filename, "w")) == NULL) {
806  std::fprintf(stderr, "Error: Cannot open file: %s\n", filename);
807  return 0;
808  }
809  } else
810  out_file = stdout;
811 
812  if (Ptrfmt != NULL) strcpy(Ptrfmt, "(8I10)");
813  ParseIfmt(Ptrfmt, &Ptrperline, &Ptrwidth);
814  std::sprintf(pformat, "%%%dd", Ptrwidth);
815  ptrcrd = (N + 1) / Ptrperline;
816  if ((N + 1) % Ptrperline != 0) ptrcrd++;
817 
818  if (Indfmt == NULL) Indfmt = Ptrfmt;
819  ParseIfmt(Indfmt, &Indperline, &Indwidth);
820  std::sprintf(iformat, "%%%dd", Indwidth);
821  indcrd = nz / Indperline;
822  if (nz % Indperline != 0) indcrd++;
823 
824  if (Type[0] != 'P') { /* Skip if pattern only */
825  if (Valfmt != NULL) strcpy(Valfmt, "(4E20.13)");
826  ParseRfmt(Valfmt, &Valperline, &Valwidth, &Valprec, &Valflag);
827  if (Valflag == 'D') *std::strchr(Valfmt, 'D') = 'E';
828  if (Valflag == 'F')
829  std::sprintf(vformat, "%% %d.%df", Valwidth, Valprec);
830  else
831  std::sprintf(vformat, "%% %d.%dE", Valwidth, Valprec);
832  valcrd = nvalentries / Valperline;
833  if (nvalentries % Valperline != 0) valcrd++;
834  } else
835  valcrd = 0;
836 
837  if (Nrhs > 0) {
838  if (Rhsfmt == NULL) Rhsfmt = Valfmt;
839  ParseRfmt(Rhsfmt, &Rhsperline, &Rhswidth, &Rhsprec, &Rhsflag);
840  if (Rhsflag == 'F')
841  std::sprintf(rformat, "%% %d.%df", Rhswidth, Rhsprec);
842  else
843  std::sprintf(rformat, "%% %d.%dE", Rhswidth, Rhsprec);
844  if (Rhsflag == 'D') *std::strchr(Rhsfmt, 'D') = 'E';
845  rhscrd = nrhsentries / Rhsperline;
846  if (nrhsentries % Rhsperline != 0) rhscrd++;
847  if (Rhstype[1] == 'G') rhscrd += rhscrd;
848  if (Rhstype[2] == 'X') rhscrd += rhscrd;
849  rhscrd *= Nrhs;
850  } else
851  rhscrd = 0;
852 
853  totcrd = 4 + ptrcrd + indcrd + valcrd + rhscrd;
854 
855  /* Print header information: */
856 
857  std::fprintf(out_file, "%-72s%-8s\n%14d%14d%14d%14d%14d\n", Title, Key, totcrd,
858  ptrcrd, indcrd, valcrd, rhscrd);
859  std::fprintf(out_file, "%3s%11s%14d%14d%14d\n", Type, " ", M, N, nz);
860  std::fprintf(out_file, "%-16s%-16s%-20s", Ptrfmt, Indfmt, Valfmt);
861  if (Nrhs != 0) {
862  /* Print Rhsfmt on fourth line and */
863  /* optional fifth header line for auxillary vector information: */
864  std::fprintf(out_file, "%-20s\n%-14s%d\n", Rhsfmt, Rhstype, Nrhs);
865  } else
866  std::fprintf(out_file, "\n");
867 
868  offset = 1 - _SP_base; /* if base 0 storage is declared (via macro definition), */
869  /* then storage entries are offset by 1 */
870 
871  /* Print column pointers: */
872  for (i = 0; i < N + 1; i++) {
873  entry = colptr[i] + offset;
874  std::fprintf(out_file, pformat, entry);
875  if ((i + 1) % Ptrperline == 0) std::fprintf(out_file, "\n");
876  }
877 
878  if ((N + 1) % Ptrperline != 0) std::fprintf(out_file, "\n");
879 
880  /* Print row indices: */
881  for (i = 0; i < nz; i++) {
882  entry = rowind[i] + offset;
883  std::fprintf(out_file, iformat, entry);
884  if ((i + 1) % Indperline == 0) std::fprintf(out_file, "\n");
885  }
886 
887  if (nz % Indperline != 0) std::fprintf(out_file, "\n");
888 
889  /* Print values: */
890 
891  if (Type[0] != 'P') { /* Skip if pattern only */
892 
893  for (i = 0; i < nvalentries; i++) {
894  std::fprintf(out_file, vformat, val[i]);
895  if ((i + 1) % Valperline == 0) std::fprintf(out_file, "\n");
896  }
897 
898  if (nvalentries % Valperline != 0) std::fprintf(out_file, "\n");
899 
900  /* If available, print right hand sides,
901  guess vectors and exact solution vectors: */
902  acount = 1;
903  linemod = 0;
904  if (Nrhs > 0) {
905  for (i = 0; i < Nrhs; i++) {
906  for (j = 0; j < nrhsentries; j++) {
907  std::fprintf(out_file, rformat, rhs[j]);
908  if (acount++ % Rhsperline == linemod) std::fprintf(out_file, "\n");
909  }
910  if ((acount - 1) % Rhsperline != linemod) {
911  std::fprintf(out_file, "\n");
912  linemod = (acount - 1) % Rhsperline;
913  }
914  rhs += nrhsentries;
915  if (Rhstype[1] == 'G') {
916  for (j = 0; j < nrhsentries; j++) {
917  std::fprintf(out_file, rformat, guess[j]);
918  if (acount++ % Rhsperline == linemod) std::fprintf(out_file, "\n");
919  }
920  if ((acount - 1) % Rhsperline != linemod) {
921  std::fprintf(out_file, "\n");
922  linemod = (acount - 1) % Rhsperline;
923  }
924  guess += nrhsentries;
925  }
926  if (Rhstype[2] == 'X') {
927  for (j = 0; j < nrhsentries; j++) {
928  std::fprintf(out_file, rformat, exact[j]);
929  if (acount++ % Rhsperline == linemod) std::fprintf(out_file, "\n");
930  }
931  if ((acount - 1) % Rhsperline != linemod) {
932  std::fprintf(out_file, "\n");
933  linemod = (acount - 1) % Rhsperline;
934  }
935  exact += nrhsentries;
936  }
937  }
938  }
939  }
940 
941  if (std::fclose(out_file) != 0) {
942  std::fprintf(stderr, "Error closing file in writeHB_mat_double().\n");
943  return 0;
944  } else
945  return 1;
946 }
947 
948 int readHB_mat_char(const char* filename, int colptr[], int rowind[],
949  char val[], char* Valfmt) {
950  /****************************************************************************/
951  /* This function opens and reads the specified file, interpreting its */
952  /* contents as a sparse matrix stored in the Harwell/Boeing standard */
953  /* format and creating compressed column storage scheme vectors to hold */
954  /* the index and nonzero value information. */
955  /* */
956  /* ---------- */
957  /* **CAVEAT** */
958  /* ---------- */
959  /* Parsing real formats from Fortran is tricky, and this file reader */
960  /* does not claim to be foolproof. It has been tested for cases when */
961  /* the real values are printed consistently and evenly spaced on each */
962  /* line, with Fixed (F), and Exponential (E or D) formats. */
963  /* */
964  /* ** If the input file does not adhere to the H/B format, the ** */
965  /* ** results will be unpredictable. ** */
966  /* */
967  /****************************************************************************/
968  std::FILE* in_file;
969  int i, j, ind, col, offset, count, last;
970  int Nrow, Ncol, Nnzero, Nentries, Nrhs;
971  int Ptrcrd, Indcrd, Valcrd, Rhscrd;
972  int Ptrperline, Ptrwidth, Indperline, Indwidth;
973  int Valperline, Valwidth, Valprec;
974  int Valflag; /* Indicates 'E','D', or 'F' float format */
975  char* ThisElement;
976  char line[BUFSIZ];
977  char Title[73], Key[9], Type[4] = "XXX", Rhstype[4];
978  char Ptrfmt[17], Indfmt[17], Rhsfmt[21];
979 
980  if ((in_file = std::fopen(filename, "r")) == NULL) {
981  std::fprintf(stderr, "Error: Cannot open file: %s\n", filename);
982  return 0;
983  }
984 
985  readHB_header(in_file, Title, Key, Type, &Nrow, &Ncol, &Nnzero, &Nrhs,
986  Ptrfmt, Indfmt, Valfmt, Rhsfmt,
987  &Ptrcrd, &Indcrd, &Valcrd, &Rhscrd, Rhstype);
988 
989  /* Parse the array input formats from Line 3 of HB file */
990  ParseIfmt(Ptrfmt, &Ptrperline, &Ptrwidth);
991  ParseIfmt(Indfmt, &Indperline, &Indwidth);
992  if (Type[0] != 'P') { /* Skip if pattern only */
993  ParseRfmt(Valfmt, &Valperline, &Valwidth, &Valprec, &Valflag);
994  if (Valflag == 'D') {
995  *std::strchr(Valfmt, 'D') = 'E';
996  }
997  }
998 
999  /* Read column pointer array: */
1000 
1001  offset = 1 - _SP_base; /* if base 0 storage is declared (via macro definition), */
1002  /* then storage entries are offset by 1 */
1003 
1004  ThisElement = (char*)malloc(Ptrwidth + 1);
1005  if (ThisElement == NULL) IOHBTerminate("Insufficient memory for ThisElement.");
1006  *(ThisElement + Ptrwidth) = '\0';
1007  count = 0;
1008  for (i = 0; i < Ptrcrd; i++) {
1009  if (std::fgets(line, BUFSIZ, in_file) == NULL) {
1010  std::fprintf(stderr, "Error: Failed to read from file.\n");
1011  return 0;
1012  }
1013  if (std::sscanf(line, "%*s") < 0)
1014  IOHBTerminate("Trilinos_Util_iohb.cpp: Null (or blank) line in pointer data region of HB file.\n");
1015  col = 0;
1016  for (ind = 0; ind < Ptrperline; ind++) {
1017  if (count > Ncol) break;
1018  std::strncpy(ThisElement, line + col, Ptrwidth);
1019  /*ThisElement = substr(line,col,Ptrwidth);*/
1020  colptr[count] = std::atoi(ThisElement) - offset;
1021  count++;
1022  col += Ptrwidth;
1023  }
1024  }
1025  free(ThisElement);
1026 
1027  /* Read row index array: */
1028 
1029  ThisElement = (char*)malloc(Indwidth + 1);
1030  if (ThisElement == NULL) IOHBTerminate("Insufficient memory for ThisElement.");
1031  *(ThisElement + Indwidth) = '\0';
1032  count = 0;
1033  for (i = 0; i < Indcrd; i++) {
1034  if (std::fgets(line, BUFSIZ, in_file) == NULL) {
1035  std::fprintf(stderr, "Error: Failed to read from file.\n");
1036  return 0;
1037  }
1038  if (std::sscanf(line, "%*s") < 0)
1039  IOHBTerminate("Trilinos_Util_iohb.cpp: Null (or blank) line in index data region of HB file.\n");
1040  col = 0;
1041  for (ind = 0; ind < Indperline; ind++) {
1042  if (count == Nnzero) break;
1043  std::strncpy(ThisElement, line + col, Indwidth);
1044  /*ThisElement = substr(line,col,Indwidth);*/
1045  rowind[count] = std::atoi(ThisElement) - offset;
1046  count++;
1047  col += Indwidth;
1048  }
1049  }
1050  free(ThisElement);
1051 
1052  /* Read array of values: AS CHARACTERS*/
1053 
1054  if (Type[0] != 'P') { /* Skip if pattern only */
1055 
1056  if (Type[0] == 'C')
1057  Nentries = 2 * Nnzero;
1058  else
1059  Nentries = Nnzero;
1060 
1061  ThisElement = (char*)malloc(Valwidth + 1);
1062  if (ThisElement == NULL) IOHBTerminate("Insufficient memory for ThisElement.");
1063  *(ThisElement + Valwidth) = '\0';
1064  count = 0;
1065  for (i = 0; i < Valcrd; i++) {
1066  if (std::fgets(line, BUFSIZ, in_file) == NULL) {
1067  std::fprintf(stderr, "Error: Failed to read from file.\n");
1068  return 0;
1069  }
1070  if (std::sscanf(line, "%*s") < 0)
1071  IOHBTerminate("Trilinos_Util_iohb.cpp: Null (or blank) line in value data region of HB file.\n");
1072  if (Valflag == 'D') {
1073  while (std::strchr(line, 'D')) *std::strchr(line, 'D') = 'E';
1074  }
1075  col = 0;
1076  for (ind = 0; ind < Valperline; ind++) {
1077  if (count == Nentries) break;
1078  ThisElement = &val[count * Valwidth];
1079  std::strncpy(ThisElement, line + col, Valwidth);
1080  /*std::strncpy(ThisElement,substr(line,col,Valwidth),Valwidth);*/
1081  if (Valflag != 'F' && std::strchr(ThisElement, 'E') == NULL) {
1082  /* insert a char prefix for exp */
1083  last = std::strlen(ThisElement);
1084  for (j = last + 1; j >= 0; j--) {
1085  ThisElement[j] = ThisElement[j - 1];
1086  if (ThisElement[j] == '+' || ThisElement[j] == '-') {
1087  ThisElement[j - 1] = Valflag;
1088  break;
1089  }
1090  }
1091  }
1092  count++;
1093  col += Valwidth;
1094  }
1095  }
1096  }
1097 
1098  return 1;
1099 }
1100 
1101 int readHB_newmat_char(const char* filename, int* M, int* N, int* nonzeros, int** colptr,
1102  int** rowind, char** val, char** Valfmt) {
1103  std::FILE* in_file;
1104  int Nrhs;
1105  int Ptrcrd, Indcrd, Valcrd, Rhscrd;
1106  int Valperline, Valwidth, Valprec;
1107  int Valflag; /* Indicates 'E','D', or 'F' float format */
1108  char Title[73], Key[9], Type[4] = "XXX", Rhstype[4];
1109  char Ptrfmt[17], Indfmt[17], Rhsfmt[21];
1110 
1111  if ((in_file = std::fopen(filename, "r")) == NULL) {
1112  std::fprintf(stderr, "Error: Cannot open file: %s\n", filename);
1113  return 0;
1114  }
1115 
1116  *Valfmt = (char*)malloc(21 * sizeof(char));
1117  if (*Valfmt == NULL) IOHBTerminate("Insufficient memory for Valfmt.");
1118  readHB_header(in_file, Title, Key, Type, M, N, nonzeros, &Nrhs,
1119  Ptrfmt, Indfmt, (*Valfmt), Rhsfmt,
1120  &Ptrcrd, &Indcrd, &Valcrd, &Rhscrd, Rhstype);
1121  std::fclose(in_file);
1122  ParseRfmt(*Valfmt, &Valperline, &Valwidth, &Valprec, &Valflag);
1123 
1124  *colptr = (int*)malloc((*N + 1) * sizeof(int));
1125  if (*colptr == NULL) IOHBTerminate("Insufficient memory for colptr.\n");
1126  *rowind = (int*)malloc(*nonzeros * sizeof(int));
1127  if (*rowind == NULL) IOHBTerminate("Insufficient memory for rowind.\n");
1128  if (Type[0] == 'C') {
1129  /*
1130  std::fprintf(stderr, "Warning: Reading complex data from HB file %s.\n",filename);
1131  std::fprintf(stderr, " Real and imaginary parts will be interlaced in val[].\n");
1132  */
1133  /* Malloc enough space for real AND imaginary parts of val[] */
1134  *val = (char*)malloc(*nonzeros * Valwidth * sizeof(char) * 2);
1135  if (*val == NULL) IOHBTerminate("Insufficient memory for val.\n");
1136  } else {
1137  if (Type[0] != 'P') {
1138  /* Malloc enough space for real array val[] */
1139  *val = (char*)malloc(*nonzeros * Valwidth * sizeof(char));
1140  if (*val == NULL) IOHBTerminate("Insufficient memory for val.\n");
1141  }
1142  } /* No val[] space needed if pattern only */
1143  return readHB_mat_char(filename, *colptr, *rowind, *val, *Valfmt);
1144 }
1145 
1146 int readHB_aux_char(const char* filename, const char AuxType, char b[]) {
1147  /****************************************************************************/
1148  /* This function opens and reads the specified file, placing auxilary */
1149  /* vector(s) of the given type (if available) in b : */
1150  /* Return value is the number of vectors successfully read. */
1151  /* */
1152  /* AuxType = 'F' full right-hand-side vector(s) */
1153  /* AuxType = 'G' initial Guess vector(s) */
1154  /* AuxType = 'X' eXact solution vector(s) */
1155  /* */
1156  /* ---------- */
1157  /* **CAVEAT** */
1158  /* ---------- */
1159  /* Parsing real formats from Fortran is tricky, and this file reader */
1160  /* does not claim to be foolproof. It has been tested for cases when */
1161  /* the real values are printed consistently and evenly spaced on each */
1162  /* line, with Fixed (F), and Exponential (E or D) formats. */
1163  /* */
1164  /* ** If the input file does not adhere to the H/B format, the ** */
1165  /* ** results will be unpredictable. ** */
1166  /* */
1167  /****************************************************************************/
1168  std::FILE* in_file;
1169  int i, j, n, maxcol, start, stride, col, last, linel, nvecs, rhsi;
1170  int Nrow, Ncol, Nnzero, Nentries, Nrhs;
1171  int Ptrcrd, Indcrd, Valcrd, Rhscrd;
1172  int Rhsperline, Rhswidth, Rhsprec;
1173  int Rhsflag;
1174  char Title[73], Key[9], Type[4] = "XXX", Rhstype[4];
1175  char Ptrfmt[17], Indfmt[17], Valfmt[21], Rhsfmt[21];
1176  char line[BUFSIZ];
1177  char* ThisElement;
1178 
1179  if ((in_file = std::fopen(filename, "r")) == NULL) {
1180  std::fprintf(stderr, "Error: Cannot open file: %s\n", filename);
1181  return 0;
1182  }
1183 
1184  readHB_header(in_file, Title, Key, Type, &Nrow, &Ncol, &Nnzero, &Nrhs,
1185  Ptrfmt, Indfmt, Valfmt, Rhsfmt,
1186  &Ptrcrd, &Indcrd, &Valcrd, &Rhscrd, Rhstype);
1187 
1188  if (Nrhs <= 0) {
1189  std::fprintf(stderr, "Warn: Attempt to read auxillary vector(s) when none are present.\n");
1190  return 0;
1191  }
1192  if (Rhstype[0] != 'F') {
1193  std::fprintf(stderr, "Warn: Attempt to read auxillary vector(s) which are not stored in Full form.\n");
1194  std::fprintf(stderr, " Rhs must be specified as full. \n");
1195  return 0;
1196  }
1197 
1198  /* If reading complex data, allow for interleaved real and imaginary values. */
1199  if (Type[0] == 'C') {
1200  Nentries = 2 * Nrow;
1201  } else {
1202  Nentries = Nrow;
1203  }
1204 
1205  nvecs = 1;
1206 
1207  if (Rhstype[1] == 'G') nvecs++;
1208  if (Rhstype[2] == 'X') nvecs++;
1209 
1210  if (AuxType == 'G' && Rhstype[1] != 'G') {
1211  std::fprintf(stderr, "Warn: Attempt to read auxillary Guess vector(s) when none are present.\n");
1212  return 0;
1213  }
1214  if (AuxType == 'X' && Rhstype[2] != 'X') {
1215  std::fprintf(stderr, "Warn: Attempt to read auxillary eXact solution vector(s) when none are present.\n");
1216  return 0;
1217  }
1218 
1219  ParseRfmt(Rhsfmt, &Rhsperline, &Rhswidth, &Rhsprec, &Rhsflag);
1220  maxcol = Rhsperline * Rhswidth;
1221 
1222  /* Lines to skip before starting to read RHS values... */
1223  n = Ptrcrd + Indcrd + Valcrd;
1224 
1225  for (i = 0; i < n; i++) {
1226  if (std::fgets(line, BUFSIZ, in_file) == NULL) {
1227  std::fprintf(stderr, "Error: Failed to read from file.\n");
1228  return 0;
1229  }
1230  }
1231 
1232  /* start - number of initial aux vector entries to skip */
1233  /* to reach first vector requested */
1234  /* stride - number of aux vector entries to skip between */
1235  /* requested vectors */
1236  if (AuxType == 'F')
1237  start = 0;
1238  else if (AuxType == 'G')
1239  start = Nentries;
1240  else
1241  start = (nvecs - 1) * Nentries;
1242  stride = (nvecs - 1) * Nentries;
1243 
1244  if (std::fgets(line, BUFSIZ, in_file) == NULL) {
1245  std::fprintf(stderr, "Error: Failed to read from file.\n");
1246  return 0;
1247  }
1248  linel = std::strchr(line, '\n') - line;
1249  if (std::sscanf(line, "%*s") < 0)
1250  IOHBTerminate("Trilinos_Util_iohb.cpp: Null (or blank) line in auxillary vector data region of HB file.\n");
1251  col = 0;
1252  /* Skip to initial offset */
1253 
1254  for (i = 0; i < start; i++) {
1255  col += Rhswidth;
1256  if (col >= (maxcol < linel ? maxcol : linel)) {
1257  if (std::fgets(line, BUFSIZ, in_file) == NULL) {
1258  std::fprintf(stderr, "Error: Failed to read from file.\n");
1259  return 0;
1260  }
1261  linel = std::strchr(line, '\n') - line;
1262  if (std::sscanf(line, "%*s") < 0)
1263  IOHBTerminate("Trilinos_Util_iohb.cpp: Null (or blank) line in auxillary vector data region of HB file.\n");
1264  col = 0;
1265  }
1266  }
1267 
1268  if (Rhsflag == 'D') {
1269  while (std::strchr(line, 'D')) *std::strchr(line, 'D') = 'E';
1270  }
1271  /* Read a vector of desired type, then skip to next */
1272  /* repeating to fill Nrhs vectors */
1273 
1274  for (rhsi = 0; rhsi < Nrhs; rhsi++) {
1275  for (i = 0; i < Nentries; i++) {
1276  if (col >= (maxcol < linel ? maxcol : linel)) {
1277  if (std::fgets(line, BUFSIZ, in_file) == NULL) {
1278  std::fprintf(stderr, "Error: Failed to read from file.\n");
1279  return 0;
1280  }
1281  linel = std::strchr(line, '\n') - line;
1282  if (std::sscanf(line, "%*s") < 0)
1283  IOHBTerminate("Trilinos_Util_iohb.cpp: Null (or blank) line in auxillary vector data region of HB file.\n");
1284  if (Rhsflag == 'D') {
1285  while (std::strchr(line, 'D')) *std::strchr(line, 'D') = 'E';
1286  }
1287  col = 0;
1288  }
1289  ThisElement = &b[i * Rhswidth];
1290  std::strncpy(ThisElement, line + col, Rhswidth);
1291  if (Rhsflag != 'F' && std::strchr(ThisElement, 'E') == NULL) {
1292  /* insert a char prefix for exp */
1293  last = std::strlen(ThisElement);
1294  for (j = last + 1; j >= 0; j--) {
1295  ThisElement[j] = ThisElement[j - 1];
1296  if (ThisElement[j] == '+' || ThisElement[j] == '-') {
1297  ThisElement[j - 1] = Rhsflag;
1298  break;
1299  }
1300  }
1301  }
1302  col += Rhswidth;
1303  }
1304  b += Nentries * Rhswidth;
1305 
1306  /* Skip any interleaved Guess/eXact vectors */
1307 
1308  for (i = 0; i < stride; i++) {
1309  col += Rhswidth;
1310  if (col >= (maxcol < linel ? maxcol : linel)) {
1311  if (std::fgets(line, BUFSIZ, in_file) == NULL) {
1312  std::fprintf(stderr, "Error: Failed to read from file.\n");
1313  return 0;
1314  }
1315  linel = std::strchr(line, '\n') - line;
1316  if (std::sscanf(line, "%*s") < 0)
1317  IOHBTerminate("Trilinos_Util_iohb.cpp: Null (or blank) line in auxillary vector data region of HB file.\n");
1318  col = 0;
1319  }
1320  }
1321  }
1322 
1323  std::fclose(in_file);
1324  return Nrhs;
1325 }
1326 
1327 int readHB_newaux_char(const char* filename, const char AuxType, char** b, char** Rhsfmt) {
1328  std::FILE* in_file;
1329  int Ptrcrd, Indcrd, Valcrd, Rhscrd;
1330  int Nrow, Ncol, Nnzero, Nrhs;
1331  int Rhsperline, Rhswidth, Rhsprec;
1332  int Rhsflag;
1333  char Title[73], Key[9], Type[4] = "XXX", Rhstype[4];
1334  char Ptrfmt[17], Indfmt[17], Valfmt[21];
1335 
1336  if ((in_file = std::fopen(filename, "r")) == NULL) {
1337  std::fprintf(stderr, "Error: Cannot open file: %s\n", filename);
1338  return 0;
1339  }
1340 
1341  *Rhsfmt = (char*)malloc(21 * sizeof(char));
1342  if (*Rhsfmt == NULL) IOHBTerminate("Insufficient memory for Rhsfmt.");
1343  readHB_header(in_file, Title, Key, Type, &Nrow, &Ncol, &Nnzero, &Nrhs,
1344  Ptrfmt, Indfmt, Valfmt, (*Rhsfmt),
1345  &Ptrcrd, &Indcrd, &Valcrd, &Rhscrd, Rhstype);
1346  std::fclose(in_file);
1347  if (Nrhs == 0) {
1348  std::fprintf(stderr, "Warn: Requested read of aux vector(s) when none are present.\n");
1349  return 0;
1350  } else {
1351  ParseRfmt(*Rhsfmt, &Rhsperline, &Rhswidth, &Rhsprec, &Rhsflag);
1352  if (Type[0] == 'C') {
1353  std::fprintf(stderr, "Warning: Reading complex aux vector(s) from HB file %s.", filename);
1354  std::fprintf(stderr, " Real and imaginary parts will be interlaced in b[].");
1355  *b = (char*)malloc(Nrow * Nrhs * Rhswidth * sizeof(char) * 2);
1356  if (*b == NULL) IOHBTerminate("Insufficient memory for rhs.\n");
1357  return readHB_aux_char(filename, AuxType, *b);
1358  } else {
1359  *b = (char*)malloc(Nrow * Nrhs * Rhswidth * sizeof(char));
1360  if (*b == NULL) IOHBTerminate("Insufficient memory for rhs.\n");
1361  return readHB_aux_char(filename, AuxType, *b);
1362  }
1363  }
1364 }
1365 
1366 int writeHB_mat_char(const char* filename, int M, int N,
1367  int nz, const int colptr[], const int rowind[],
1368  const char val[], int Nrhs, const char rhs[],
1369  const char guess[], const char exact[],
1370  const char* Title, const char* Key, const char* Type,
1371  char* Ptrfmt, char* Indfmt, char* Valfmt, char* Rhsfmt,
1372  const char* Rhstype) {
1373  /****************************************************************************/
1374  /* The writeHB function opens the named file and writes the specified */
1375  /* matrix and optional right-hand-side(s) to that file in Harwell-Boeing */
1376  /* format. */
1377  /* */
1378  /* For a description of the Harwell Boeing standard, see: */
1379  /* Duff, et al., ACM TOMS Vol.15, No.1, March 1989 */
1380  /* */
1381  /****************************************************************************/
1382  std::FILE* out_file;
1383  int i, j, acount, linemod, entry, offset;
1384  int totcrd, ptrcrd, indcrd, valcrd, rhscrd;
1385  int nvalentries, nrhsentries;
1386  int Ptrperline, Ptrwidth, Indperline, Indwidth;
1387  int Rhsperline, Rhswidth, Rhsprec;
1388  int Rhsflag;
1389  int Valperline, Valwidth, Valprec;
1390  int Valflag; /* Indicates 'E','D', or 'F' float format */
1391  char pformat[16], iformat[16], vformat[19], rformat[19];
1392 
1393  if (Type[0] == 'C') {
1394  nvalentries = 2 * nz;
1395  nrhsentries = 2 * M;
1396  } else {
1397  nvalentries = nz;
1398  nrhsentries = M;
1399  }
1400 
1401  if (filename != NULL) {
1402  if ((out_file = std::fopen(filename, "w")) == NULL) {
1403  std::fprintf(stderr, "Error: Cannot open file: %s\n", filename);
1404  return 0;
1405  }
1406  } else
1407  out_file = stdout;
1408 
1409  if (Ptrfmt != NULL) strcpy(Ptrfmt, "(8I10)");
1410  ParseIfmt(Ptrfmt, &Ptrperline, &Ptrwidth);
1411  std::sprintf(pformat, "%%%dd", Ptrwidth);
1412 
1413  if (Indfmt == NULL) Indfmt = Ptrfmt;
1414  ParseIfmt(Indfmt, &Indperline, &Indwidth);
1415  std::sprintf(iformat, "%%%dd", Indwidth);
1416 
1417  if (Type[0] != 'P') { /* Skip if pattern only */
1418  if (Valfmt != NULL) strcpy(Valfmt, "(4E20.13)");
1419  ParseRfmt(Valfmt, &Valperline, &Valwidth, &Valprec, &Valflag);
1420  std::sprintf(vformat, "%%%ds", Valwidth);
1421  }
1422 
1423  ptrcrd = (N + 1) / Ptrperline;
1424  if ((N + 1) % Ptrperline != 0) ptrcrd++;
1425 
1426  indcrd = nz / Indperline;
1427  if (nz % Indperline != 0) indcrd++;
1428 
1429  valcrd = nvalentries / Valperline;
1430  if (nvalentries % Valperline != 0) valcrd++;
1431 
1432  if (Nrhs > 0) {
1433  if (Rhsfmt == NULL) Rhsfmt = Valfmt;
1434  ParseRfmt(Rhsfmt, &Rhsperline, &Rhswidth, &Rhsprec, &Rhsflag);
1435  std::sprintf(rformat, "%%%ds", Rhswidth);
1436  rhscrd = nrhsentries / Rhsperline;
1437  if (nrhsentries % Rhsperline != 0) rhscrd++;
1438  if (Rhstype[1] == 'G') rhscrd += rhscrd;
1439  if (Rhstype[2] == 'X') rhscrd += rhscrd;
1440  rhscrd *= Nrhs;
1441  } else
1442  rhscrd = 0;
1443 
1444  totcrd = 4 + ptrcrd + indcrd + valcrd + rhscrd;
1445 
1446  /* Print header information: */
1447 
1448  std::fprintf(out_file, "%-72s%-8s\n%14d%14d%14d%14d%14d\n", Title, Key, totcrd,
1449  ptrcrd, indcrd, valcrd, rhscrd);
1450  std::fprintf(out_file, "%3s%11s%14d%14d%14d\n", Type, " ", M, N, nz);
1451  std::fprintf(out_file, "%-16s%-16s%-20s", Ptrfmt, Indfmt, Valfmt);
1452  if (Nrhs != 0) {
1453  /* Print Rhsfmt on fourth line and */
1454  /* optional fifth header line for auxillary vector information: */
1455  std::fprintf(out_file, "%-20s\n%-14s%d\n", Rhsfmt, Rhstype, Nrhs);
1456  } else
1457  std::fprintf(out_file, "\n");
1458 
1459  offset = 1 - _SP_base; /* if base 0 storage is declared (via macro definition), */
1460  /* then storage entries are offset by 1 */
1461 
1462  /* Print column pointers: */
1463  for (i = 0; i < N + 1; i++) {
1464  entry = colptr[i] + offset;
1465  std::fprintf(out_file, pformat, entry);
1466  if ((i + 1) % Ptrperline == 0) std::fprintf(out_file, "\n");
1467  }
1468 
1469  if ((N + 1) % Ptrperline != 0) std::fprintf(out_file, "\n");
1470 
1471  /* Print row indices: */
1472  for (i = 0; i < nz; i++) {
1473  entry = rowind[i] + offset;
1474  std::fprintf(out_file, iformat, entry);
1475  if ((i + 1) % Indperline == 0) std::fprintf(out_file, "\n");
1476  }
1477 
1478  if (nz % Indperline != 0) std::fprintf(out_file, "\n");
1479 
1480  /* Print values: */
1481 
1482  if (Type[0] != 'P') { /* Skip if pattern only */
1483  for (i = 0; i < nvalentries; i++) {
1484  std::fprintf(out_file, vformat, val + i * Valwidth);
1485  if ((i + 1) % Valperline == 0) std::fprintf(out_file, "\n");
1486  }
1487 
1488  if (nvalentries % Valperline != 0) std::fprintf(out_file, "\n");
1489 
1490  /* Print right hand sides: */
1491  acount = 1;
1492  linemod = 0;
1493  if (Nrhs > 0) {
1494  for (j = 0; j < Nrhs; j++) {
1495  for (i = 0; i < nrhsentries; i++) {
1496  std::fprintf(out_file, rformat, rhs + i * Rhswidth);
1497  if (acount++ % Rhsperline == linemod) std::fprintf(out_file, "\n");
1498  }
1499  if (acount % Rhsperline != linemod) {
1500  std::fprintf(out_file, "\n");
1501  linemod = (acount - 1) % Rhsperline;
1502  }
1503  if (Rhstype[1] == 'G') {
1504  for (i = 0; i < nrhsentries; i++) {
1505  std::fprintf(out_file, rformat, guess + i * Rhswidth);
1506  if (acount++ % Rhsperline == linemod) std::fprintf(out_file, "\n");
1507  }
1508  if (acount % Rhsperline != linemod) {
1509  std::fprintf(out_file, "\n");
1510  linemod = (acount - 1) % Rhsperline;
1511  }
1512  }
1513  if (Rhstype[2] == 'X') {
1514  for (i = 0; i < nrhsentries; i++) {
1515  std::fprintf(out_file, rformat, exact + i * Rhswidth);
1516  if (acount++ % Rhsperline == linemod) std::fprintf(out_file, "\n");
1517  }
1518  if (acount % Rhsperline != linemod) {
1519  std::fprintf(out_file, "\n");
1520  linemod = (acount - 1) % Rhsperline;
1521  }
1522  }
1523  }
1524  }
1525  }
1526 
1527  if (std::fclose(out_file) != 0) {
1528  std::fprintf(stderr, "Error closing file in writeHB_mat_char().\n");
1529  return 0;
1530  } else
1531  return 1;
1532 }
1533 
1534 int ParseIfmt(char* fmt, int* perline, int* width) {
1535  /*************************************************/
1536  /* Parse an *integer* format field to determine */
1537  /* width and number of elements per line. */
1538  /*************************************************/
1539  char* tmp;
1540  if (fmt == NULL) {
1541  *perline = 0;
1542  *width = 0;
1543  return 0;
1544  }
1545  upcase(fmt);
1546  tmp = std::strchr(fmt, '(');
1547  tmp = substr(fmt, tmp - fmt + 1, std::strchr(fmt, 'I') - tmp - 1);
1548  *perline = std::atoi(tmp);
1549  if (*perline == 0) *perline = 1;
1550  if (tmp != NULL) free((void*)tmp);
1551  tmp = std::strchr(fmt, 'I');
1552  tmp = substr(fmt, tmp - fmt + 1, std::strchr(fmt, ')') - tmp - 1);
1553  *width = std::atoi(tmp);
1554  if (tmp != NULL) free((void*)tmp);
1555  return *width;
1556 }
1557 
1558 int ParseRfmt(char* fmt, int* perline, int* width, int* prec, int* flag) {
1559  /*************************************************/
1560  /* Parse a *real* format field to determine */
1561  /* width and number of elements per line. */
1562  /* Also sets flag indicating 'E' 'F' 'P' or 'D' */
1563  /* format. */
1564  /*************************************************/
1565  char* tmp;
1566  char* tmp1;
1567  char* tmp2;
1568  char* tmp3;
1569  int len;
1570 
1571  if (fmt == NULL) {
1572  *perline = 0;
1573  *width = 0;
1574  flag = NULL;
1575  return 0;
1576  }
1577 
1578  upcase(fmt);
1579  if (std::strchr(fmt, '(') != NULL) fmt = std::strchr(fmt, '(');
1580  if (std::strchr(fmt, ')') != NULL) {
1581  tmp2 = std::strchr(fmt, ')');
1582  while (std::strchr(tmp2 + 1, ')') != NULL) {
1583  tmp2 = std::strchr(tmp2 + 1, ')');
1584  }
1585  *(tmp2 + 1) = '\0';
1586  }
1587  if (std::strchr(fmt, 'P') != NULL) /* Remove any scaling factor, which */
1588  { /* affects output only, not input */
1589  if (std::strchr(fmt, '(') != NULL) {
1590  tmp = std::strchr(fmt, 'P');
1591  if (*(++tmp) == ',') tmp++;
1592  tmp3 = std::strchr(fmt, '(') + 1;
1593  len = tmp - tmp3;
1594  tmp2 = tmp3;
1595  while (*(tmp2 + len) != '\0') {
1596  *tmp2 = *(tmp2 + len);
1597  tmp2++;
1598  }
1599  *(std::strchr(fmt, ')') + 1) = '\0';
1600  }
1601  }
1602  if (std::strchr(fmt, 'E') != NULL) {
1603  *flag = 'E';
1604  } else if (std::strchr(fmt, 'D') != NULL) {
1605  *flag = 'D';
1606  } else if (std::strchr(fmt, 'F') != NULL) {
1607  *flag = 'F';
1608  } else {
1609  std::fprintf(stderr, "Real format %s in H/B file not supported.\n", fmt);
1610  return 0;
1611  }
1612  tmp = std::strchr(fmt, '(');
1613  tmp = substr(fmt, tmp - fmt + 1, std::strchr(fmt, *flag) - tmp - 1);
1614  *perline = std::atoi(tmp);
1615  if (*perline == 0) *perline = 1;
1616  if (tmp != NULL) free((void*)tmp);
1617  tmp = std::strchr(fmt, *flag);
1618  if (std::strchr(fmt, '.')) {
1619  tmp1 = substr(fmt, std::strchr(fmt, '.') - fmt + 1, std::strchr(fmt, ')') - std::strchr(fmt, '.') - 1);
1620  *prec = std::atoi(tmp1);
1621  if (tmp1 != NULL) free((void*)tmp1);
1622  tmp1 = substr(fmt, tmp - fmt + 1, std::strchr(fmt, '.') - tmp - 1);
1623  } else {
1624  tmp1 = substr(fmt, tmp - fmt + 1, std::strchr(fmt, ')') - tmp - 1);
1625  }
1626  *width = std::atoi(tmp1);
1627  if (tmp1 != NULL) free((void*)tmp1);
1628  return *width;
1629 }
1630 
1631 char* substr(const char* S, const int pos, const int len) {
1632  int i;
1633  char* SubS;
1634  if ((size_t)pos + len <= std::strlen(S)) {
1635  SubS = (char*)malloc(len + 1);
1636  if (SubS == NULL) IOHBTerminate("Insufficient memory for SubS.");
1637  for (i = 0; i < len; i++) SubS[i] = S[pos + i];
1638  SubS[len] = '\0';
1639  } else {
1640  SubS = NULL;
1641  }
1642  return SubS;
1643 }
1644 
1645 #include <cctype>
1646 void upcase(char* S) {
1647  /* Convert S to uppercase */
1648  int i, len;
1649  len = std::strlen(S);
1650  for (i = 0; i < len; i++)
1651  S[i] = std::toupper(S[i]);
1652 }
1653 
1654 void IOHBTerminate(const char* message) {
1655  std::fprintf(stderr, "%s", message);
1656  std::exit(1);
1657 }
void start()
Start the deep_copy counter.